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ABSTFL\CT 

An  acoustic  ray  tracing  algorithm  is  developed  and  coupled  with  a 
thermodynamic  upper  ocean  mixed  layer  model.  For  a  test  case,  the  coupled  niixed 
layer-acoustic  model  is  applied  to  a  specific  area  in  the  western  Mediterranean  Sea. 
Climatological  atmospheric  forcing  is  used  to  provide  boundary  conditions  for  the 
mixed  layer  for  short  periods  of  time  (from  few  hours  to  three  days)  during  different 
seasons.  The  response  of  the  acoustic  model  to  the  predicted  changes  in  the  sound- 
speed  profile  is  analyzed  to  show  dependence  of  acoustic  propagation  upon  the  surface 
atmospheric  forcing  and  the  season.  The  atmospheric  factors  such  as  wind,  rain,  and 
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factor  which  influences  the  acoustic  propagation.  The  effect  of  heavy  rain  with  light 
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I.  INTRODUCTION 

An  acoustic  ray  tracing  program  is  coupled  with  the  Oceanic  Boundan.'  Layer 
Model  developed  by  Garwood  [Ref.  1].  The  OBL  model  is  a  one-dimensional,  second- 
order  turbulence  closure,  vertically  integrated  model  of  the  ocean  surface  turbulent 
boundary  layer,  usmg  a  two-component  turbulent  kinetic  energy  budget  with  a  mean 
turbulent  field  closure. 

Fisher  [Ref  2]  investigated  the  variabihty  and  sensitivity  of  a  coupled  model 
system.  He  found  that  the  OBL  model,  when  integrated  in  time  at  a  single  point 
(Ocean  Station  Papa  50° N  145"\V),  predicted  inixed-layer  structure  better  than  did  the 
Expanded  Ocean  Thermal  Structure  (LOTS)  system  which  was  currently  in  use  at  the 
Fleet  Numerical  Oceanography  Center  (FNOC).  McManus  [Ref.  3]  evaluated  the 
acoustic  performance  of  a  coupled  model  system  at  a  line  of  stations  in  the  northeast 
Pacific  Ocean.  In  both  cases,  the  thermodynamic  model  was  initialized  with  observed 
temperature  profiles,  and  the  surface  boundary  conditions  were  given  by  the  currently 
available  meteorological  informations.  Then,  the  thermodynamic  forecasts  were  input 
into  acoustic  models,  such  as  RAYMODE  or  FACT,  and  the  acoustic  performance 
was  analyzed  using  the  median  detection  range  (.VI DR)  and  the  convergence  zone 
range  (CZR). 

This  research  is  the  first  attempt  to  link  in  a  single  program  the  OBL  model  with 
an  acoustic  model.  SimpHcity  (compared  to  the  operational  models  available  in  the 
U.S.  Navy)  and  classification  restrictions  lead  us  to  develop  in  Chapter  II  an  algorithm 
for  acoustic  ray  tracing.  As  no  such  routine  vv^as  available  at  NFS,  a  copy  is  attached 
in  the  appendix  for  further  use  by  students  of  the  Air-Ocean  Sciences  Department.  This 
simple  subroutine  allows  the  iniluence  of  the  atmospheric  forcing  on  the  underwater 
sound  propagation  to  be  qualitatively  analyzed.  A  summary  of  the  leading  principles 
and  equations  of  the  OBL  model  is  given  in  Chapter  III.  Chapter  IV  develops  m  detail 
the  actual  coupUng  of  the  two  models  into  a  single  computer  code.  Chapter  V  gives  an 
example  of  the  use  of  this  coupled  model  applied  to  a  specific  area  in  the  western 
.Mediterranean  Sea,  for  dilTerent  periods  of  the  year  having  particular  acoustic 
properties.  All  the  simulations  were  integrated  in  time  using  climatological  data  over 
short  time  periods  varying  from  a  few  hours  to  three  days. 
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II.  ACOUSTIC  MODEL 

A.       INTRODUCTION 

Since  we  assumed  in  this  research  thai  the  ocean  is  horizontally  stratified,  the 
temperature  T  and  the  salinity  S  are  only  functions  of  the  depth  }■  and  cannot  vary 
with  the  range  z.  Hence,  the  speed  of  sound  c  is  only  a  function  of  depth  y.  In  that 
case,  according  to  Ziomek  [Ref  4:  p. 236],  the  general  form  of  the  equation  for  the 
horizontal  range  travelled  by  an  acoustic  ray  is  the  following  : 

2  =  ^0  +  ^I''   TTT^v^'A  (2-1) 

where  b  is  the  ray  parameter  and  is  given  by 

b  =  sinPCyg)  .'  c(yo),  (2.2) 

(y^.z^)  are  the  coordinates  of  the  source,  and  Plyg)  is  the  initial  angle  of  propagation. 

Thus,  theoretically,  knowing  the  sound-speed  profile  versus  the  depth,  we  can 
plot  the  curve  giving  the  path  of  an  acoustic  ray. 

We  chose  not  to  use  equation  (2.1)  for  the  following  reasons: 

According  to  Snell's  law. 


sinpfy)   _  smP(y,,)_  ^  ^^ 

^(y)  ciy^) 

and.  at  a  turning  point.  PCy^p)  =  ^'2  and  c(y.  )  =  1/b  such  that  the  denominator  of 
the  integrand  of  (2.1)  goes  to  zero  and  the  integration  cannot  be  carried  out.  Thus 
(2.1)  is  only  vahd  between  two  turning  points. 

Numerically,  the  integration  of  (2.1)  would  be  carried  out  between  all  the  turning 
points  of  the  ray  by  discretizing  the  sound-speed  profile  c(y)  with  as  small  a  depth 
increment  as  required  to  get  an  acceptable  result.  Furthermore,  during  the  integration, 
we  would  have  to  test  for  the  occurence  of  a  turning  point. 
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However,  a  sound-speed  profile  can  be  approximated  by  straight  line  segments 
matching  the  profile  as  best  as  possible.  The  more  segments  chosen,  the  better  the 
computation.  Roughly,  we  often  choose  10  to  15  segments  depending  on  the  shape  of 
the  profile.  Then,  as  explained  in  the  next  section,  the  integration  of  (2.1)  can  be 
solved  analytically  for  each  segment  of  straight  line  and  leads  to  very  simple  equations. 


0                        2.                                           Z 

z 

V 

y. 

Po 

^ 

""^-^-^ 

p 

y> 

r 

Figure  2.1     Ray  path  confined  to  the  YZ  plane. 

We  shall  keep  in  mind  that  this  acoustic  model  has  to  be  coupled  with  an 
Oceanic  Planetary  Boundary  Layer  (OPBL)  model  [Ref  1]  which  resolves  the 
temperature,  salinity,  and  depth  of  the  upper  mixed  layer  of  the  ocean  using  a  one- 
meter-step  increment  for  the  depth. 

These  are  the  reasons  why  we  chose  to  discretize  the  sound-speed  profile  with  a 
one-meter-depth  increment,  assume  the  profile  to  be  a  straight  line  segment  within  each 
depth  increment,  and  we  use  the  equations  derived  in  the  next  section  instead  of 
numerically  evaluating  the  integral  of  (2.1).  We  found  this  method  easier  for  handling 
the  turning  point  problem. 

B.   SPEED  OF  SOUND  AS  A  FUNCTION  OF  DEPTH  WITH  CONSTANT 
GRADIENT 

The  sound-speed  profile  is  given  by 


c(y)  =  c.  +  g  X  (y-y.) 


(2.4) 
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where  g  is  a  constant  (with  units  ofs    ''  referred  to  as  the  gradient  since 


defy) 


=  8- 


(2.5) 


Starting  with  Snell's  law  (2.3),  differentiation  leads  to: 


co.sp(y)  dp  =  b  dc(y)  =  bg  dy. 


(2.6) 


Referring  to  Figure  2.2,  it  can  be  easily  seen  that 


dv 

■^  =  cosp(y). 


(2.7) 


p 

dz 

z 

>. 

^(y) 

\as 

^ 

ay 

Vy 

f 

Figure  2.2     An  infinitesimal  element  of  arc  length  ds 
at  an  arbitrary  point  P  along  a  ray  path  in  YZ  plane. 


Substituting  (2.7)  in  (2.6)  yields 


dP(y) 

ds 


=  bg. 


(2.8) 
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Equation  (2.S)  indicates  that  the  curvature  along  a  ray  path  is  constant.    Thus,  the  ray 
path  is  an  arc  of  a  circle.    From  Figure  2.2,  we  have 


dz  =  tanP(y)  dy  =  sinP(y)  -^— —  .  (2.9) 

cosply) 

Using  (2.6),  we  get 

dz  =  sinp(v)'^^  .  (2.10) 

'      bg 

Integration  leads  to 

where  p^  =  P(yg)  is  the  initial  angle  of  propagation  at  the  source,  and  P  =  p(y)  is  the 
angle  of  propagation  at  position  P.   This  leads  to 


z  =  Zq  +  -^(cospQ-cosP).  (2.12) 

Solving  for  the  depth  y  in  (2.4)  yields 


y  =  Yo  +  — (c(y)-  Cq)  ,  (2.13) 


and.  using  Snell's  law  (2.3),  we  get 

y  =  y^  +  l_(sinp-sinpo).  (2.14) 

We  can  verify  that  (2.12)  and  (2.14)  are  the  parametric  equations  (parameter  p)  of  the 
circle  siven  bv  : 
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y  -  (Vq  -  T-^  -  +  z  +  (zg  +  -^)  -  =  (— r  (2.15) 

u     bg  ^       bg  bg 

centered  at  (  y^  -  sinp^  bg  .  -{z^  -r  cosP^  bg)  )  and  having  a  radius  of  1  bg. 

C.       .ALGORITHM 

For  this  research,  we  used  a  modiiled  version  of  the  basic  [-""ORTR.AN  subroutine 
R.'\.Y  given  in  appendix  A.  In  this  section  we  are  going  to  analyze  in  detail  the 
different  parts  of  the  algorithm  used  in  the  subroutine  Ri\Y. 

Tlie  goal  o[  this  subroutine  is  to  plot  the  acoustic  ray  path,  that  is,  the  curve 
giving  the  range  z  in  kilometers  versus  the  depth  y  in  meters  of  an  acoustic  ray 
emanating  from  a  source  at  a  depth  y^  at  a  given  initial  angle  of  propagation.  The  ray 
path  plot  is  presented  beside  the  graph  of  the  sound-speed  profile  c{y)  in  m  sec  versus 
depth  y  in  meters. 

As  we  mentioned  previously,  we  used  a  one-meter-step  depth  increment  and  kept 
track  of  the  depth  all  along  the  ray  by  using  an  increasing  or  decreasing  index  k 
depending  on  the  ray  going  upward  or  downward.  The  index  k  is,  in  fact,  the  integer 
value  of  the  depth. 

The  angles  of  propagation  P  and  Pq  are  referenced  from  the  Y  axis.  For  plotting 
purposes,  we  defined  the  depth  vector  {  y-,  i  =  0,NN  }  and  the  range  vector  (  z,  i  = 
O.NN  }.   At  the  i^^  point  of  the  ray  path,  the  depth  index  k  is 

k  =  y^,  (2.16) 

the  "initial"  angle  oC  propagation  is  Pg,  the  sound-speed  is  Cn  and,  its  gradient,  gu,  is 
2iven  bv 


gk  =  '^k+l-^k-  ('-1') 

At  the  (i-t-1)^"  point  on  the  ray  path,  the  depth  index  is  k±I  the  "final"  angle  of 
propagation  is  p,  the  sound-speed  is  cj^^  j  and  the  gradient  is  gj^^  j  .  An  illustration 
of  these  different  parameters  is  given  in  Figure  2.3. 
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For  that  case,  according  to  (2.3)  and  (2.12),  we  would  have 


}\+l  =  k-1 


(2.18) 


and 


=    Z:    + 


^1+1 "  -^ib^: 


(cosPp,  -  cosP) 


k-1 


(2.19) 


with 


sinP   _  sinP^. 


''k-1        ^source 


nurce  =  ^  . 


(2.20) 


0  C.  Ck     Ck.  C 

0       1^, 

2;„ 

z 

k.l 

1 
1 

^ 

e. 

/ 

\ 

^u. 

y" 

^^ 

K 

\ 

^^L<--''^ 

"" 

\ 

yc 

-^p. 

\ 

y> 

r 

' 

Figure  2.3     Parameters  describing  a  one-meter-depth  increment 

positive  gradient  sound-speed  profile 

for  an  upward  travelling  acoustic  ray  path. 
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Now  we  are  going  to  examine  in  detail  the  dillerent  cases  occuring  in  the 
subroutine  leading  to  some  IF. ..THEN. ..ELSE...  statements,  depending  on  the  sign  of 
the  gradient  (or  the  curvature)  and  the  upward  or  downward  direction  of  the  ray  at  the 
i     point. 

1.  Case  I  :  g,..|>0,  ?i^>0. 

At  the  i^     point  of  the  ray  path  P:,  corresponding  to  the  depth  index  k,  the 
gradient  above  is  (see  Figure  2.4) 


8k-l  =  ^k-^k-1 


(2.21) 


and  the  gradient  below  is 


^k  =  ^k+l-<^k 


(2.22) 


Figure  2.4    Case  1  :  gj^.|  >  0,  gj^>0. 

We  have  three  difTercnt  cases  to  consider  depending  on  the  value  of  P^  with 
regard  to  7r/2  and  P^.  The  critical  angle  p^,  is  the  value  of  the  initial  angle  of 
propagation  at  P-  leading  to  a  final  angle  P  =  n;2  one  meter  deeper,  that  is,  leaduig  to 
a  turning  point  one  meter  deeper. 
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According  to  Sncll's  law, 


Ci, 


sinp^  =  —^  (2.23; 


C  r 


or 


P^  =  arcsin  ^^-  .  (2.24) 

^k+1 

Let  us  examine  the  three  cases  Pq  >  71/2,  Pq  <  P^  and  P^  <  Pq  <  ti/2  since  a  strict 
equality  does  not  apply  for  real  numbers  in  FORTRAN  language. 

a.    pQ  >   7T/2  : 

According  to  (2,3)  and  since  the  ARCSIN  function  gives  a  result  between  0 
and  71/2,  the  angle  of  propagation  at  P-^_  |  is  : 

P  =  7C  -  arcsin(bc^_|)  .  (2.25) 

Substituing  for  P  in  (2.12),  we  get  the  coordinates  of  the  next  point  ?[+  [  : 

yi+l  =  k-l  (2.26) 

and 

z-^_i  =  z-  +  - —   (cospQ  +  cos[arcsin(bC|^.|)])  (2.27) 

with 

b  =  s^^Psource  '  Source  •  (2.28) 
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*•    P„  <  P,  : 

This  time,  using  the  same  equations  (2.3)  and  (2.12)  at  the  point  ^[+1-  '^^'c 
have  : 

P  =  arcsin(bci,  +  i),  (2.29) 


Yi+l  =  k+1  ,  (2.30) 

and 

^1+  1  "  ^1  ^  "gr"  ^^°^Po  "  cos[arcsin(bcj^+  |)]}  .  (2.31) 

In  this  case,  we  have  to  deal  with  a  turning  point  within  the  segment 

Let  us  derive  the  general  formula  for  computing  z-_)_  ^  when  a  turning  point 
is  encountered.  This  formula  will  be  used  when  we  study  the  cases  generated  by  the 
other  sign  possibilities  of  the  gradient  g. 

Applying  (2.12)  first  from  Pq  to  P',  then  from  P'  to  P  we  have  (see  Figure 
2.5)  : 


z'  =  Zn  +  — — -  (cosPn  -  cos7r/2)  (2.32) 


and 


z'  +  [COSTT,  2  -  cos(7i  -  Pn)]  ,  (2.33) 

b2 


or 


19 


z    =  z. 


bg 


cosp. 


(2.34) 


and 


z  =  z'  -I ^^  cospQ 


]_ 

bg 


(2.35) 


Substituing  (2.34)  in  (2.35)  yields 


^  ^  ^0  ^  "bT  '°'^0  • 


(2.36) 


Figure  2.5    Turning  point  treatment. 


of  P 


Going  back  to  our  case  of  interest  and  using  (2.36),  we  get  the  coordinates 


i+l 


■14-1  =  Vi  =  k 


(2.37) 
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and 


^i+l  =  ^i  + 


bgj. 


cosp. 


(2.38) 


2.  Case2:g,..j<0,  gj.<0 

This  case  is  ver>'  similar  to  case  1,  but  now  the  critical  angle  p^  is  given  by  : 


P^  =  ;:  -  arcsin(cj^  C|^_j) 


(2.39) 


Following  the  notation  of  Figure  2.6  and  applying  the  same  basic  equations  (2.3), 
(2.12),  and  (2.36),  we  have  three  new  cases  to  examine. 


CVvi     Ck   Ck.,    C         Zl 


k.t 


y 


yi 


c<^.<^c 


z 

— > 


'j.^.^U^Vr    P^^, 


Figure  2.6    Case  2  :  g^_|  <0,  gj^<0. 


a,   %  <  Ti'l  : 


P  =  arcsin(bcj^^_  j)  , 


(2.40) 
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Yi+I  =  k+1,  (2.41) 

and 

Zj+  1  =  Zj  +  -^  [cosPq  -  cos[arcsin{bcj^^  j)]}  .  (2.42) 

b.  p,  >  8,  : 

P  =  7T  -  arcsin(bcj^.j)  ,  (2.43) 


y,^l  =  k-[,  (2.44) 

and 

Z;j_i  =  z-  + (cosp,.  +  cos[arcsin(bC].  i)])  .  (2.45) 

C.     TT  2    <    Po    <    Pc  : 

y.^^  =  y-  =  k  .  (2.46) 

and 

^1+1  ==^1^  lJ-,^°^Po-  (2.47) 

»k-l 

3.  Case  3:  gk.i>0,  gi^<0 

According  to  Figure  2.7,  we  now  have  only  two  separate  cases  to  consider, 
depending  on  the  value  of  Pq  with  regard  to  k:2.  Returning  to  the  two  previous  cases, 
we  will  find  some  similarities. 


apply. 


apply. 


a.  p(5  >  71/2  : 

This  case  is  similar  to  case   l.a,  and  equations  (2.25),  (2.26),  and  (2.27) 

b.  Po  <  71/2  : 

This  case  is  similar  to  case  2. a,  and  equations  (2.40),  (2.41),  and  (2.42) 


Ck-. 

Ck,.  Ck    c 

— > 

1                1     1   ^ 

;           1    ; 
1              '    I 

-  -• 

P  .  ^ ...      u  . 

\ 

\    '     '^'-' 

^^-^ 

f 

/      ^' 

>P^ 

y> 

' 

X             Ji.** 

Figure  2.7    Case  3  :  gj^_  j  >  0,  gj^  <  0. 

4.  Case  4:  gk.i<0,  gk>0 

According  to  Figure  2.8,  we  now  have  four  difTerent  cases  to  examine  with 
two  possibilities  for  P^  . 


Pj,  =  7t  -  arcsin 


A. 


Ci. 


(2.48) 


k-1 


or 


arcsm 


C|. 


(2.49) 


k+1 
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The  four  different  cases  shown  above  reduce  to  cases  previously  analyzed  in 
sections  for  case  1  and  case  2. 


Figure  2.8    Case  4  :  gj^. ^  <  0,  gj^ >  0. 

a-    Po  >   Pc  ^ 

This  case  is  similar  to  case  2.b,  and  equations  (2.43),  (2.44),  and  (2.45) 


apply. 


b.   P,  >  Po  >  nil  : 

This  case  is  similar  to  case  2.c,  and  equations  (2.46)  and  (2.47)  apply. 

This  case  is  similar  to  case  l.b,  and  equations  (2.29),  (2.30),  and  (2.31) 
apply. 

d.  nil  >  Po  >  p,.  : 

This  case  corresponds  to  case  l.c,  and  equations  (2.37)  and  (2.38)  apply. 
This  last  case  ends  the  cascade  of  IF. ..THEN. ..ELSE...  statements  that  we  used  to 
write  this  subroutine. 
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5.  Surface  and  bottom  reflections 

The  last  point  we  need  to  discuss  is  how  to  handle  surface  and  bottom 
reflections,  that  is,  when  the  depth  index  k.  reaches  the  values  of  0  or  N.  The  value  N 
corresponds  to  the  maximum  depth  given  in  the  profile  and  is  assumed  to  be  the  depth 
of  the  sea  floor.  In  the  case  where  the  sound-speed  profile  does  not  extend  to  the 
bottom,  it  is  always  possible  to  extrapolate  the  sound-speed  profile  to  the  bottom  using 
the  gradient  of  the  profile  corresponding  to  the  last  straight  line  segment  or,  using  an 
average  value  given  by  the  climatology. 


l,_o 

p^.:>\ 

P... 

/ 

9o   -^ 

■^ 

§^5f- — — - 

^ 

% 

f^-.i 

^v^ 

PCH 

^^ 

f- 

See    cases     1.1  ar^^  1.5 

i^-.o 

ft  ^V 

9o 

r\ 

k-.i 

/         \ 

R„ 

5<?e     case  2.1 

\ 

Figure  2.9     Perfect  surface  reflection  possibilities. 

a.   Surface  reflection 

In  this  simple  acoustic  model,  we  assume  a  perfect  surface  reflection.  Thus, 
when  the  depth  index  k  reaches  the  value  0,  we  just  have  to  maintain  symmetry  (i.e., 
the  angle  of  reflection  is  equal  to  the  angle  of  incidence)  with  regard  to  the  horizontal 
to  get  the  correct  angle  Pq  before  applying  the  same  equations  as  in  cases  l.b,  I.c  and 
2. a  of  Fieure  2.9  . 
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b.   Bottom  reflection 

Following  the  same  idea,  we  assume  a  perfect  bottom  reflection.  Thus 
when  the  depth  index  k  reaches  its  maximum  value  N,  we  maintain  symmetry'  (i.e.,  the 
angle  of  reflection  is  equal  to  the  angle  of  incidence)  with  regard  to  the  horizontal  (see 
Figure  2.10). 

Notice  that  the  fictitious  gradients  g_j  and  %-^  have  been  created  in  the 
program  to  handle  the  cascade  of  IF. ..THEN. ..ELSE...  statements. 


brH.I        - 
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jN.t 

\ 

X             /^ 
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\                                       / 
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fe-.N   _ 

See 
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cases 
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Figure  2.10     Perfect  bottom  reflection  possibilities. 


D.       SUMMARY  AND  DIRECTIVES  TO  USE  THE  SUBROUTINE  RAY 

1.  Main  characteristics 

The  computation  in  the  subroutine  RAY  is  based  on  representing  an  arbitrary 
sound-speed  profile  by  straight  line  segments  in  one-meter  depth  increments  and  thus 
allows  the  use  of  any  complicated  sound-speed  profile  to  obtain  a  precise  and  smooth 
acoustic  ray  trace. 
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A  plot  of  the  sound-speed  profile  is  provided  at  the  same  time  as  the  plot  of 
the  bundle  of  rays. 

The  main  assumptions  underlying  this  program  arc  the  following  : 

•  perfect  surface  reflection, 

•  perfect  bottoni  reflection, 

•  flat  bottom. 

2.   Features  of  this  subroutine 

This  subroutine  uses  many  DISSPL.-X  graphics  statements.  Hence,  the  main 
program  calling  subroutine  I^\Y  has  to  be  executed  using  the  command  DISSPLA.  A 
choice  of  multiple  display  devices  is  provided  to  the  user  through  the  use  of  comment 
cards  in  front  of  the  statem.ents  CALL  COMPRS,  CALL  TEK618,  or  CALL 
CX-41(4I07). 

The  variables  required  when  using  R.'W  are  the  following  : 

YO  :  depth  of  the  source  in  meters. 

M  :  number  of  ray  path  desired, 

BET  :  array  (size  M)  of  the  initial  angles  of  propagation  P<;Qurrf;  in  degrees, 

BO  :  initial  angle  of  the  upper  ray  of  the  bundle, 

DB  :  increment  of  initial  angles  in  BET, 

M\H- 1  :  number  of  points  in  the  provided  sound-speed  profile. 

CC  :  provided  values  of  sound  speed  in  m.sec, 

YY  :  corresponding  depth  in  meters, 

R-ANGE  :  maximum  range  desired  in  kilometers, 

XN  :  index  of  range, 

Y  :  array  (size  XN)  of  depth  in  meters, 

Z  :  array  (size  XX)  of  horizontal  distances  in  kilometers, 

X  :  integer  value  of  the  depth  of  the  sea  floor  in  meters. 

C  :  array  (size  X)  of  sound  speed  (m  sec), 

G  :  array  (size  X)  of  its  gradient  (sec    ), 

YC  :  array  (size  X)  of  depths  (meters). 
The  bundle  of  rays  to  be  traced  is  defined  by  the  number  of  rays  M,  the  initial 
angle  BO  of  the  upper  ray  of  the  bundle  ,  and  the  angle  increment  DB  between  two 
rays.  Inside  the  subroutine,  we  can  also  set  different  equations  to  define  a  bundle  of 
rays.  And  finally,  by  deleting  these  few  lines  involved  with  these  com.putations,  we  can 
provide  our  own  array  BET  of  initial  angles  PsQ^J■ce■ 


No  relation  has  been  derived  between  NX  and  the  maximum  range  FLXNGE 
since  this  relation  depends  on  many  ditTercnt  parameters  such  as  the  Pjource'^'  ^^'^^ 
sound-speed  profile,  and  the  maximum  range.  If  during  a  plot  a  ray  ends  before 
reaching  the  right  side  of  the  graph  (especially  the  steepest  rays),  a  higher  number  NN 
has  to  be  set.  For  instance,  a  value  NN=  5000  seems  to  handle  most  of  the  reasonable 
steep  PjQuj-ce  S^^'^^  '^'^^  '^  maximum  range  of  30  km. 

The  arrays  Y  and  Z  have  been  created  for  computing  and  plotting  the  ray 
paths  and  the  arrays  C,  G  and  YC  are  for  plotting  the  sound-speed  profile. 
3.  Example 

A  result  of  the  R.-XY  subroutine  is  displayed  in  Figure  2.11.  The  provided 
parameters  and  sound-speed  profile  are  listed  in  Table  1  . 


TABLE  1 

PAFL'\METERS  AND  DATA  USED  IN  THE  EXAMPLE 

OF  SUBROUTINE  RAY 

SHOWN  IN  FIGURE2.il 

^'0 

30. 

YY 

0. 

cc 

1482. 

M 

10 

30. 

1482.9 

NN 

5000 

55. 

1484.1 

R.ANGE 

30. 

90. 

1480.5 

N 

300 

150. 

1482.0 

MM 

6 

225. 
300. 

1484.1 
1486.2 
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(W)    Hld3a 


Figure  2.11     Ray  tracing  provided  by  subroutine  RAY 
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III.  OCEAN  MIXED  LAYER  MODEL 

A.       INTRODUCTION 

The  study  of  the  oceanic  turbulent  boundan'  layer  is  a  relatively  recent  field  in 
Physical  Oceanography.  The  model  used  for  this  research  has  been  developed  by 
Garwood  [Refs.  1,5].  These  papers  give  the  formulation  of  a  unified  mathematical 
model  of  the  one-dimensional,  non-stationary,  oceanic  turbulent  boundary  layer. 

The  study  of  these  top  few  tens  of  meters  of  the  ocean  is  of  considerable 
scientific  interest.  It  influences  and  can  be  related  to  the  general  circulation  of  the 
ocean  [Ref  6].  The  thermal  structure  associated  with  this  boundary  layer  must  be 
considered  when  making  medium  and  long-range  weather  forecast,  since  a  large  part  of 
the  atmospheric  energy  supply  comes  from  the  heat  exchanged  with  the  ocean.  This 
layer  is  also  a  region  of  primar>'  biological  productivity,  which  is  of  significant 
ecological  and  economic  importance.  Finally,  as  an  important  military  application,  this 
study  can  be  used  in  the  modehng  of  acoustic  propagation  in  the  ocean,  which  is  the 
goal  of  this  research. 

1 .  Characteristics  of  the  oceanic  mixed  layer 

The  oceanic  mixed  layer  is  that  fully  turbulent  region  of  the  upper  ocean  that 
is  bounded  above  by  the  air-sea  interface  and  below  by  a  dynamicaly  stable  water 
mass.  The  wind  and  intermittent  upward  surface  buoyancy  flux  through  the  surface 
(surface  cooling  at  night  or  in  winter  time  for  example)  are  the  sources  of  mechanical 
energy  for  the  generation  of  these  turbulences.  Figure  3.1  gives  a  general  picture 
idealizing  density  and  mean  velocity  profiles  of  the  ocean  mixed  layer. 

2.  Generalities  on  the  dynamics  of  the  mixing 

The  depth  of  the  ocean  wind-mixed  surface  layer  is  typically  on  the  order  of 
ten  to  one  hundred  meters.  The  horizontal  scale  size  is  that  of  the  internal  Rossby 
radius,  typically  20  to  50  km.  These  two  dominant  scale  sizes  are  much  smaller  than 
the  horizontal  scale  size  of  the  driving  meteorological  disturbances,  water  mass 
features,  and  distances  to  lateral  boundaries.    The  approximation  of  local  horizontal 


In  this  chapter,  in  order  to  be  consistent  with  the  notation  used  in  [Refs.  1.5]  as 
in  most  of  the  geophysical  sciences  publications,  we  chose  z  to  be  the  upward  vertical 
axis  and  (x,y)  to  be  the  horizontal  coordinates  as  depicted  in  Figure  3.1.  The  vertical 
component  of  the  fluid  velocity  will  be  w  and  the  horizontal  components  v^-ill  be  u  and 
v. 
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homogeneity  for  all  mean  variables  is  usually  accurate  and  is  a  basic  assumption  in  this 
model. 
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Figure  3.1     Idealized  model  for  ocean  mixed  layer.. 

A  sharp  density  discontinuity  of  thickness  5  (see  Figure  3.1)  separates  the 
layer  from  a  stable  non  turbulent  thermocline.  Minimal  stress  at  the  bottom  together 
with  high  turbulence  intensity  leads  also  to  an  approximate  vertical  uniformity  in  mean 
velocity  and  density.  We  shall  note  that  only  small  gradients  in  these  mean  variables 
give  rise  to  large  turbulent  fluxes. 

The  mechanical  energy  budget  for  the  ocean  mixed  layer  is  depicted  in  Figure 
3.2.  Deepening  of  the  mixed  layer  is  accomplished  by  entrainment  of  the  more  dense 
underlying  water  into  the  turbulent  region  above.  This  process  leads  to  a  potential 
energy  increase  and  cannot  take  place  without  an  energy  source:  the  turbulent  kinetic 
energy  of  the  mixed  layer  above. 

Retreat  occurs  when  the  vertical  component  of  the  turbulence  is  insufficient  to 
transport  heat,  momentum,  and  turbulence  to  an  earlier-established  depth  of  mixing. 
This  happens  when  downward  heat  flux  (buoyant  damping)  and  dissipation  effects 
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Figure  3.2     Mechanical  energy  budget  for  the  ocean  mixed  layer. 
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exceed  the  wind-stirring  shear-production  of  turbulence,  creating  a  stratification  of  the 
upper  ocean. 

Thermal  energy  and  mechanical  energy  received  from  the  atmosphere  not  only 
control  the  local  dynamics,  but  the  layer  itself  modulates  the  flux  of  this  energy  to  the 
deeper  water  masses.  Entrainmcnt  also  converts  some  of  the  mean  flow  energy  into 
turbulent  energy,  over  and  above  the  wmd-stress  production. 

Finally,  substantial  barotropic  and  baroclinic  features,  such  as  tidal  motion 
and  internal  waves,  can  be  linearly  superimposed.    The  mean  fields  of  concern  arc 
therefore  the  horizontally  homogeneous  components  of  the  total  fields. 
3.  The  model  and  its  features  not  previously  demonstrated 

The  vertical  and  horizontal  components  of  the  turbulent  kinetic  energy  (TKE) 
are  determined  implicitely.  along  with  layer  depth,  mean  momentum,  and  mean 
buoyancy.    Layer  growth  and  retreat  are  predicted. 

Specific  features  of  the  TKE  budget  include  mean  turbulent  field  modehng  of 
the  dissipation  term,  the  energy  redistribution  term,  and  the  term  for  the  convergence 
of  buoyancy  fiux  at  the  stable  interface  as  shown  in  a  following  section.  Then  an 
entrainmcnt  hypothesis  dependent  upon  the  relative  distribution  of  the  TKE  between 
horizontal  and  vertical  components  permits  the  closure  of  the  system  of  equations. 

The  m.odel  differs  from  earUer  models  in  the  following  ways.  First,  the 
amount  of  wind  generated  TKE  to  be  used  in  mixing  is  a  function  of  the  ratio  of  the 
mixed  layer  depth  to  the  Obukhov  mixing  length  L.  Second,  viscous  dissipation  is 
dependent  on  a  local  Rossby  number.  Finally,  separate  vertical  and  horizontal 
equations  for  TKE  arc  used,  permitting  a  more  consistent  interpretation  of  both 
entraining  and  retreating  mixed  layers. 

B.        FORMULATION  OF  THE  EQUATIONS  USED  IN  THIS  OBL  MODEL 

1.  Generalities 

The  underlying  principles  employed  in  studying  the  mixed  layer  are  the 
combined  conservation  of  mass,  momentum,  thermal  energy,  and  mechanical  energy. 
Most  of  the  physics  behind  such  a  one-dimensional  model  is  based  on  the  fiux  form  of 
the  Navier-Stokes  equations  with  the  Boussinesq  approximation.  The  horizontal 
homogeneity  mentioned  previously  permits  the  neglect  of  the  horizontal  gradients  of 
the  mean  fields. 

Conservation  of  buoyancy  is  employed  as  a  generalization  of  the  conservation 
of  heat  alone.  The  buoyancy  equation  is  generated  from  the  heat  and  salt  equations 
together  with  an  equation  of  state, 
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p  =  p,  [1  -a(G-eo)  +  P(s-Sq)].  (3.1) 

and  the  definition  for  buoyancy, 

b  =  glPy-p)  ;  Pq  ,  (3.2) 

where  p^  is  a  representative  density  at  the  time  and  location  of  consideration,  g  is 
acceleration  due  to  gravity,  a  and  p  are  the  expansion  coelTicients  for  heat  and  salt. 
All  variables  are  separated  into  mean  and  fluctuating  components  : 


0 

temperature 

0  =  T  +  T' 

e 

salmity 

s  =  S  +  s' 

9 

pressure 

P  =  P  -f  p' 

• 

velocities 

u  =  U„  +  U  +  u' 

V  =  V^"  +  V  -f  V' 
\y  =  \f  +  w' 

•      buoyancy  b  =  B  +  b' 

Subscript  g  denotes  gecstrophic  components. 

2.  Mean  buoyancy  and  momentum  equations 

The  first  law  of  thermodynamics  for  an  incompressible  fluid  and  the 
conservation  of  salt  mass,  neglecting  molecular  fluxes,  lead  to  the  mean  buoyancy 
equation  : 


^B,^t  =  -dh'v^'ldz  +  agQ  ,'  p^C    .  (3.3) 

By  invoking  the  Boussinesq  approximation,  dropping  the  negligible  viscous  terms, 
assuming  incompressibility,  and  subtracting  the  geostrophic  equations  from  the  total 
momentum  equations,  we  obtain  the  mean  momentum  equations  : 

dL.di  =  t\  -  du^ldz  (3.4) 


d\';di  =  -iU  -  dv'w'Idz  .  (3.5) 
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Integration  of  (3.3),  (3.4)  and  (3.5)  over  the  entrainment  zone  from  z  =  -h-6  to  z=-h 
leads  to  the  jump  conditions  for  the  turbulent  fluxes  at  the  bottom  of  the  mixed  layer  : 


-b'\v'(-h)  =  AB  dhidi  (3.6) 


-u'\v-(-h)  =  AU  c^h  ^t  (3.7) 


-v'vv'(-h)  =  \V  dhldt  (3.8) 

wliere 

b'  =  agT'  -  Pgs'  (3.9) 

and 

AB  =  agAT-  pgAS   .  (3.10) 

The  assumption  of  vertical  homogeneity  in  the  mixed  layer  permits  the  integration  of 
(3.3),  (3.4)  and  (3.5)  from  z  =  -h-6  to  z  =  0,  including  the  efTects  of  entrainment  stresses 
(3.7)  and  (3.8),  and  entrainment  buoyancy  flux  (3.6).  This  yields  the  bulk  relationships 
for  mean  buovancv  : 


h  d  <  B>;^t  +  AAB  ^h,^t  =  agQQ  /  PqC    -  b'vv'(O)  (3.11) 


ha<U>/at  +  AAU^h^t  =  -fh<V>  -  u'w'(O)  (3.12) 


ha<V>  5t  +  AAV  ^h  at  =  -fh<U>  -  v'w'(O)  (3.13) 
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where 


%  =  ^s-  (Qe  ^  Qh  +Qb)   i^  W'rn^  .  (3.14) 

The  value  of  Qq  is  the  net  solar  irradiance  at  the  surface,  minus  the  long  wave  back 
radiation,  minus  the  sensible  heat  flux,  minus  the  latent  heat  flux. 
The  function  A  is  the  Heaviside  step  function  : 

A  =  1     if    ch:ci  >  0  , 

A  =  0     if    dh.  Ci  <  0  . 
The  brackets  <   >  denotes  a  vertical  mean  through  the  mixed  layer, 

1      .0 

h-r  0  "'-h-0 

and       ^denotes  horizontal  mean, 

n  =  7:^,(    )dxdy.  (3.16) 

As  the  time  step  used  in  the  model  is  one  hour,  the  surface  boundar}'  conditions  are 
prescribed  hourly  : 


-u'w'(O)  =  t^(t)Po  (3.17) 


.v'w'(O)  =  yt).Po  (3.18) 


-bw'(O)  =  g  [  ps'w'(0,t)  -  aw'T'(0,t)  ]  (3.19) 

with 


T  =  p^C.U^o^  (3.20) 
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where  p.^  is  the  air  density,  C^  a  drag  coefTicient  and  U.^  the  wind  speed  at  10  meters, 
and  with 


s'\v'(0)  =  S(E  -  P)  (3.21) 

where  E-P  =  evaporation  minus  precipitation  in  cm'sec. 
3.  Mean  turbulent  kinetic  energy  equation 

Subtracting  the  scalar  product  or(u,v.w)  with  the  mean  momentum  equations 
from  the  mean  equation  for  the  total  mechanical  energy  yields  the  mean  TKE 
equation; 


1  cE  , ^U       —  dV       d         p'   E 

=  -  [u  w^  +  vw'— ]  4-  b'w'  -  —  w'(^^-)]  -  c  -  0  (3.22) 

2  CI  cz  oz  cz         pg   2 

where  (I)  (II)  (III)        (IV) 

E  =  u'-  ^  v'-  -f  w'^  .  (3.23) 

The  budgets  for  the  individual  components  of  TKE  can  also  be  formed  : 


1  dvi'-       dV      d  w'u'-         p't5u'     c        ^  ^   

__  =  .u'w'- — (-- —    ^  L  a-  ^^u'v'  -  n,u'w'  3.24 

2  ci  dz     dz      2  Pq^x      3-^2  ^        ^ 


-T-  =  -\w'- —{-— — )  +  1  — -  .  -  -  n,u  v'  (3.25) 

2  ^t  dz     cz      2     ^       Pn^v      3        ^  ' 


'3 


l^w-     — — ;     d    \v  ^    wp^       pew     c       _^  — :  .-.  ^^x 

-^-  =  b  w'  -  ^-  -—  ^  -^)  -^  -T--T  +  ^-,u  w'  3.26) 

2  ci  cz  2         Pq         PqCz      j 

(V)  (IV) 

where  fi.  =  (f2,.n-,.n^)  is  the  rotation  vector  for  the  earth. 

The  time  rate  of  change  of  TKE  is  usually  smaller  than  the  other  terms  and 
may  be  neglected.  The  term  (I)  represents  the  rate  of  mechanical  production  and  is  the 
dominant  source  of  TKE.    It  is  the  rate  of  conversion  o[  mean  to  TKE  bv  the  dov/n- 
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gradient  turbulent  flux  of  momentum.  The  term  (11)  represents  the  buoyancy  flux  and 
can  be  either  a  source  or  a  sink.  This  term  can  become  an  important  source  as  in  the 
case  of  strong  convective  cooling  in  the  autumn.  The  term  (III)  is  the  divergence  of 
the  turbulent  flux  of  TKE  or  turbulent  diffusion  of  itself'.  Locally,  at  the  bottom  of  the 
layer  during  occasions  of  entrainment.  a  net  convergence  of  flux  of  energy  is  necessary 
to  maintain  the  downward  buoyancy  flux  for  deepening  the  mixed  layer.  The  term  (IV) 
represents  viscous  dissipation.  Because  local  isotropy  is  assumed  for  the  dissipation 
range,  e  is  divided  equally  among  the  component  budgets.  The  terms  (V)  are  ver\' 
important  terms  which  sum  to  zero  by  continuity  in  the  TKE  equation.  They  cause  a 
redistribution  of  energy  among  u'",  v''^  and  w'-.  The  terms  (VI)  are  also  redistribution 
terms,  but  they  are  due  to  rotation  of  the  earth.  In  this  study,  we  will  neglect  them 
because  of  the  usually  short  integral  time  scale  in  comparison  with  the  time  periods  of 
interest. 

Assumption  of  vertical  homogeneity  in  the  mixed  layer  permits  the  vertical 
integration  from  z  =  -h-6  to  z  =  0  of  (3.24) +  (3.25)  and  (3.26).  By  this  step  we 
introduce  a  new  variable  h.  the  depth  of  the  mixed  layer. 

Up  to  that  point,  if  we  set  q",  the  horizontal  component  of  TKE,  such  as 
q-  -  u'-  4-  v'2  ,  (3.27) 


we  have  six  variables  h,  <B>,  <U>,  <V>,  <w'->  and  <q'>  for  five  equations 
(3.11),  (3.12),  (3.13),  vertical  integral  of  (3.24)  +  (3.25)  and  (3.26).  Therefore  a  sixth 
equation  is  needed  to  close  the  system.  Besides,  we  also  need  a  suitable  modeling  of 
the  different  terms  of  the  integrated  TKE  equation.  Following  Garwood's  arguments 
[Refs.  1,5].  a  synthesis  of  these  derivations  is  given  in  the  next  sections. 
4.  Modeling  of  the  different  terms  of  the  integrated  TKE  equation 
a.   Shear  production  and  turbulent  diffusion 

The  vertical  integral  of  terms  (I)  and  (III)  may  be  combined  to  give  the  net 
"wind-generated"  rate  of  production  : 


G  =  J^  .[(I)  +  (III)]dz  (3.2S) 

-h-O 


or 


38 


G  =  m.y-  +  1,2[(AU)2  -^  (AV)-]5lv^t  (3.29 

where  ir'-  is  the  friction  velocitv  defined  as  : 


u-  =  [u'v'(O)-  +  v'w'(O)-]'-  .  (3.30) 

b.   Net  buoyant  damping 

The  integral  of  (11)  over  the  mixed  layer  gives  the  net  buoyant  damping  for 
the  whole  laver  : 


B  =  Ih.5^.I^)dz  0-31) 


or 


'0    p 


which  can  be  rewritten  as 


B  =   1,  2  h  b'w'(-h)  -  1/2  h  u*b*  (3.33) 

where  u"b*  is  the  downward  surface  flux  of  buoyancy  : 

b*u-  =  -bV(0)  +  ag  PoCp  RQJI  +  e"^  (1-^2  yh) .  ^/yh]  .  (3.34) 

The  radiation  absorption  Q(z)  has  been  modeled  as 

Q  =  yRQ^e^z  (3.35) 

where  j  is  the  extinction  coefficient  for  net  solar  radiation,  and  RQ^  is  the  short  wave 
fraction  of  the  net  solar  irradiance.  This  model  assumes  that  the  net  long-w^ave  solar 
radiation  is  absorbed  at  the  surface,  which  leads  to  : 


-b'w'(O)  =  ag  p,3Cp   [(1-R)Q3  -  (Q,  +  Qb  +  Qh)!  ^  PgS(P-E)  (3.36) 


J 
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where  {1-R)Q    corresponds  to  the  short-wave  incoming  radiation.    On  the  other  hand, 
the  short-wave  radiation  penetrates  below  the  surface  where  it  is  absorbed  following 
the  exponential  decay  (3.35). 
c.    Viscous  dissipation 

For  a  fully  turbulent  mixed  layer,  viscous  dissipation  of  the  turbulence 
occurs  primarily  in  the  small  eddies  which  are  locally  isotropic.  The  net  rate  of 
dissipation  is  modelled  as  follows: 


D  =   f       cdz  =  mj<E>-''"  +  m5th<E>  (3.37) 


5h-6 


or 


D  =  m^<E>^'-(l  +  Rq"    m3'mj  u-;<E>  ''')  (3.38) 

where  Rg  =  u*,  hf  is  a  Rossby  number  for  the  turbulent  boundary  layer.  The  first 
term  in  D  arises  from  the  fact  that  the  time  scale  of  the  largest  eddies  can  be 
proportional  to  the  mixed  layer  depth  divided  by  the  rms  turbulent  velocity  <E>  ^ 
The  second  one  comes  from  the  assumption  that,  in  deeper  boundary  layers,  planetan.' 
rotation  (time  scale  1,  f)  turns  the  mean  shear  direction  with  depth  and  thus  influences 
the  geometrical  aspects  of  the  integral  scale. 
d.    Redistribution  of  TKE 

The  vertical  integral  of  the  pressure  redistribution  term 


,.0      p'5u. 

is  an  important  source  or  sink,  term  for  the  individual  TKE  budgets,  even  though 
Rj  -^  R.,  +  R3  =  0.   The  bulk  formulation  is 

R;  =  m2<E>  '^M<E>  -  3<u^>).  (3.40) 

The   rotational   redistribution   terms    are   assumed   to    be   of  higher   order,    and   are 
neslected  in  this  studv. 
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5.  Closure  hypothesis 

Garwood  [Ref.  1]  achieves  closure  of  the  problem  by  formulating  the  following 

mtniinmcnt  equation  : 

ah.'^t  =  m_^<\v'2>  '''-<E>   '   hAB  .  (3.41) 


C.       SUMMARY  OF  MODELED  EQUATIONS 

We  have  six  variables  h,  <  B>  ,  <U>  ,  <  V>  ,  <\v'->  ,  <q->  and  a  (inal  set  of 
six  equations  : 

•  mean  momentum  equations  : 


ha<U>,^t  4-  AAU  ^h  ^t  =  -lh<V>  -  u'\v'(0)  (3.42) 


hd<V>  di  +  AAV  ^h  ^t  =  -fh<  U>  -  v'\v'(0)  (3.43) 

•  mean  buoyancy  equation  : 


h  d  <  B>  'at  +  AAB  dh;di  =  agQQ  /  PqC     -  b'\v(0)  (3.44) 

•  horizontal  integrated  TKE  component  : 

__(h<q.>)  =  m3U--3  +  ^^-^  (3.4^) 

-m2(<E>  -  3<w'^>)<E>  "' -  2mj;3   (<£>'''  +  m3,m^fh)<E> 

♦  vertical  integrated  TKE  component  : 

^-r-(h<w'^>)  =  l/2hbV(-h)-  12  h  u''^b-  +  m,(<E>  -  3  <  w'2>  )<  E> ''''        (3.46) 

2  ci  '  ^ 


m 


^  3  (<E>   ^  +  m^  m^  fh)<E: 
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•  entrainment  equation  : 

^li^t  =  m_^<w^>  -<E>  /   hAB  .  (3.47) 


The  different  constants  used  in  the  model  are  the  following  : 

p  Representative  density  of  air 

Pq         Representative  density  of  sea-water 

a.p        Coefficients  of  thermal  expansion 

g  Gravitational  acceleration 

C  Specific  heat  at  constant  pressure 

Cj         Drag  coefficient 

f  Coriolis  parameter  depending  on  latitude 

We  can  tune  the  m.odel  by  specifying  the  different  parameters  m^  m^,  m_^.  m^  as  well 
as  the  extinction  coefficient  for  net  solar  radiation  Y,  and  the  short-wave  fraction  of  the 
solar  radiation  R. 

Initial  profiles  of  temperature  T(z)  ("C  vs  cm)  and  salinity  S(z)  (g/kg  vs  cm)  have 
to  be  provided.  The  depth  increment  in  the  model  is  100  cm.  Finally,  as  the  time  step 
is  one  hour,  the  following  boundarv'  conditions  have  to  be  prescribed  hourly  : 

•  T  Wind  stress  in  dynes 'cm^ 

•  Q  Incoming  solar  radiation  in  W,'m 

•  Q^        Back  radiation  in  W/m*" 


e 


• 


Q  Latent  heat  flux 


^e 


Q,  Sensible  heat  flux 


'h 


•  E  Evaporation  in  mm/hour 

•  P  Precipitation  in  mm/hour 


The  values  of  u*  and  u*b*  can  be  computed  using  equations  (3.17).  (3.18),  (3.30), 
(3.37)  and  (3.36). 
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IV.  COL  PLING  THE  TWO  .MODELS 

A.  GENERALITIES 

The  goal  of  this  research  is  to  couple  the  Oceanic  Boundary  Layer  (OBL)  model 
derived  by  Garwood  in  1977  with  an  acoustic  ray  tracing  program  allowing  us  to 
analyze  some  effects  of  the  atmospheric  factors  on  the  oceanic  environment  and, 
therefore,  acoustic  propagation  m  the  ocean. 

The  two  previous  chapters  gave  us  a  theoretical  approach  concerning  these  two 
models.    Now,  we  are  going  to  highhght  the  coupling  of  these  two  models,  the  inputs 
that  the  coupled  model  needs,  and  the  outputs  v/e  can  get  from  it.   The  block  diagram 
of  Figure  -i.l  depicts  the  over  plan,  as  followed  in  the  next  sections. 

Initially,  we  enter  the  OBL  model  with  a  set  of  boundary  conditions,  with  some 
initial  conditions,  and  with  a  set  of  parameters  including  the  time-step  of  integration  of 
the  OBL  model  and  the  time  interval  between  each  resulting  ray  trace.  From  the 
predicted  output  of  temperature  and  salinity  profiles,  we  compute  a  sound-velocity 
profile  which  allows  us  to  trace  the  paths  of  acoustic  rays  according  to  a  set  of 
geometrical  parameters  such  as  the  source  depth,  the  maximum  range,  and  the  initial 
angle  of  a  ray." 

B.  INITIAL  CONDITIONS 

We  have  to  initialize  the  coupled  model  with  a  temperature  profile,  that  is, 
temperature  T  in  Celsius  versus  depth  z  with  a  one-meter  increment  and  a  salinity 
profile,  that  is,  salinity  S  in  g  kg  versus  depth  z  with  a  one-meter  increment.  Tliese 
profiles  can  be  obtained  from  climatological  data  [Ref  7],  as  the  ones  we  used  to 
analyze  some  applications  of  this  coupled  model,  or  from  an  XBT  (expendable 
bathythermograph)  and  a  salinity  profile  from  climatology.  We  can  also  assume  a 
constant  salinity  profile  based  on  an  average  salinity  in  the  area  studied  because  of  the 
small  effect  of  salinity  variations  on  sound-velocity  computation. 


"In  all  of  the  previous  chapters,  we  wanted  to  be  consistent  with  the  coordinate 
systems  used  in  each  reference,  that  is,  [Ref  4]  for  Chapter  II  and  [Refs.  1.5J  for 
Chapter  III.  In  Chapter  II,  y  is  the  downward  vertical  axis  and,  in  Chapter  III.  z  is 
the  upward  vertical  axis.  Since  the  coupling  is  done  through  the  output  T  and  S 
profiles  of  the  OBL  model  feeding  the  input  profiles  of  the  acoustic  model,  this 
notation  conflict  does  not  aflcct  the  coupling  by  itself  nor  the  codes  of  the  coupled 
programs. 
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Figure  4.1     Coupling  the  OBL  model  with  the  ray  tracing  program. 
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C.       BOUNDARY  CONDITIONS 

The  boundary'  conditions  or  atmospheric  forcing  factors  we  have  to  provide  are 
listed  helow. 

The  hourly  solar  radiation  flux  Q  has  to  be  given  in  W/nr  versus  time  in  hours. 
The  climatology  gives  us  generally  a  daily  average  Q...,„.  in  our  study  we  will  simulate 
a  diurnal  cycle  for  Q,  based  on  the  following  formula  : 


Q    =Q         sin(7rt/12)     if    0<t<  12  hours, 


Q3  =  o 


if    12<t<  24  hours. 


(4T) 


as  depicted  in  Figure  4.2. 


Figure  4.2     Sim.ulated  diurnal  cycle. 
A  straight-forward  integration  leads  to  the  following  relation  : 


Q         =  nO 


(4.2) 


Simulation  of  clouds  is  done  by  allowing  only  a  part  of  this  Q^  as  input  solar  radiation 
(20  to  50° 0  of  Qj  for  example)  at  the  sea  surface. 
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The  values  of  the  latent  flux  of  heat  Q^,  the  sensible  flux  of  heat  Q^  and  the  back 
radiation  Q.  have  to  be  provided  on  an  hourly  basis.  In  the  simulations  of  the  next 
chapter,  we  will  pick  some  average  values  of  the  sum  Q^  +  Qu  +  Qu  from  the 
climatology,  depending  on  what  time  of  the  year  we  are  working  with. 

The  evaporation  rate  EV  in  mm  hour  is  also  evaluated  from  climatological  data, 
and  precipitation  rate  PR  in  mm,  hour  can  be  sim.ulated.  Hea\w  rain  cases  will  be 
considered  also  because  of  some  interesting  effects. 


TABLE  2 

WIND  SPEED  AND  WIND  STRESS  CORRESPONDENCE 

Wind 

Spe 

■ed 

Wind  Stress 

knots 

m/s 

g/ cm;  sec" 

5 

2.5 

0.106 

10 

5 

0.423 

15 

7.5 

0.959 

20 

10 

1.694 

25 

12.5 

2.647 

30 

15 

3.811 

35 

17.5 

5.188 

40 

20 

6.776 

The  wind  stress  T  in  g,'cm,sec^  has  to  be  provided  hourly  also.  In  our 
simulations,  we  will  consider  both  constant  strong  wind  as  well  as  light  wind 
conditions.    Based  on  equation  (4.3) 

^  =    Pa  Cd  ^'lo'  (^-^^ 

3 
where  an  air  density  p    of  1.21  kg/m-^  and  a  drag  coefficient  C^  of  1.4x10"    have  been 

chosen.   The  correspondence  between  wind  speed  and  wind  stress  is  given  in  Table  2. 
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D.  OTHER  INPUT  PARAMETERS 

The  other  parameters  \vc  have  to  input  in  the  coupled  model  are  the  latitude  of 
the  area  studied,  the  number  of  days  NDAYS  during  which  we  want  to  integrate  the 
model,  the  mterval  of  integration  used  by  the  model  (IDIFF=  Ihour  in  this  study),  and 
the  frequency  TI  in  hours  with  which  we  want  to  output  some  graphic  results  from  the 
acoustic  ray  tracing  program. 

E.  OUTPUTS  OF  THE  OBL  MODEL 

From  the  OBL  model  we  can  get  the  hourly  values  of  the  predicted  depth, 
temperature  and  salinity  of  the  mixed  layer,  and  we  are  able  to  follow  the  evolution  of 
these  different  parameters  on  tabulated  outputs. 

The  link  with  the  acoustic  model  is  done  through  the  hourly  computed  profiles 
T(z)  and  S(z).  Depending  on  the  frequency  TI,  these  profiles  (renamed  T(y)  and  S(y) 
for  notational  consistency  as  mentioned  in  a  previous  note)  are  used  to  compute  the 
sound-speed  profile  c(y)  by  a  special  subroutine  SVEL. 

F.  SOUND-SPEED  COMPUTATION 

In  the  subroutine  computing  the  sound-speed,  we  used  the  equation  derived  by 
Chen  and  Millero  [Ref  8].  The  sound-speed  obtained  agrees  with  the  data  of  Del 
Grosso  at  one  atmosphere  [Ref  9]  and  with  the  data  of  Wilson  at  high  pressures 
[Ref  10].  Other  formulas  for  sound-speed  are  mentioned  by  Urick  [Ref  11],  but  the 
one  we  used  is  most  often  suggested  for  oceanographic  calculations.  Further  more, 
comparisons  of  graphic  outputs  do  not  reveal  any  significant  differences  bctv/een  the 
results  obtained  by  using  these  various  equations. 

G.  RAY  TRACING 

The  computed  sound-speed  profile  c(y)  is  used  as  a  main  input  in  the  acoustic 
model.  Depending  on  the  case  we  want  to  study,  a  set  of  geometrical  parameters  has 
to  be  provided  such  as  the  source  depth  YO,  the  maximum  range  R.\NGE,  and  tlie 
ma.ximum  depth  X.  The  bundle  of  rays  to  be  plotted  is  defined  by  the  number  of  rays 
M.  the  initial  angle  BO  of  the  upper  ray  of  the  bundle,  and  the  angle  increment  DEI.B 
between  two  rays.  Different  equations  for  defining  a  bundle  of  rays  can  be  set  inside 
the  subroutine.  Thus,  it  is  always  possible  to  pick  up  a  specific  bundle  of  rays  to 
highlight  some  efifects  of  the  atmospheric  forcing  on  some  precise  type  of  rays  and 
demonstrate  some  acoustic  enerev  redistribution. 
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Since,  in  this  research,  we  want  to  give  an  example  of  the  application  of  the  use 
of  this  coupled  model,  we  will  choose  generally  a  ray-bundle  width  of  5  to  8  degrees 
which  should  correspond  to  an  active  sonar  operating  at  a  frequency  in  the  range  of  5 
kHz. 

Finally,  the  main  output  of  the  coupled  model  is  a  plot  of  the  sound-speed  profile 
aside  an  acoustic  ray  tracing  for  a  given  set  of  parameters  defined  previously. 

H.      NOTE  ON  THE  ALGORITHM  USED  IN  THE  ACOUSTIC  MODEL 

Usually,  in  most  acoustic  models,  the  input  sound-speed  profile  is  approximated 
by  a  piecewisc  linear  curve  as  depicted  in  Figure  4.3. 
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Figure  4.3     Piecewise  linear  approximation 
of  a  typical  sound-speed  profile. 

Thus,  the  profile  can  be  modeled  mathematically  by  the  following  equation 


c(y)  =  S-1  "^  gi(y-yi-i) 


>'i-i  ^  y  =^  Vi 


(4.4) 
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and  is  composed  of  several  distinct  straight  line  segments,  each  with  its  own  constant 
sradicnt  c-.  In  each  such  delimited  slab  i,  because  of  the  constant  cradient  2..  the  rav 
path  is  an  arc  ot^  a  circle  for  which  an  equation  is  easy  to  compute  and  fast  to  plot. 
Since  a  profile  is  usually  divided  uito  a  few  discrete  slabs  (for  example.  10),  a  more 
computationally  elficient  subroutine  could  be  derived  from  the  one  discussed  in 
Chapter  II.  This  subroutine  would  invoke  gradient  sign  checking  and  initial  angle 
comparison  to  the  horizontal  only  when  a  ray  crosses  a  boundan.-  between  two  slabs. 

In  our  study,  this  checking  and  comparison  are  done  at  each  depth  increment  (1 
m)  as  well  as  the  plotting  of  a  new  point  belonging  to  the  ray  path.  It  is  certainly  a 
much  slower  algorithm  than  the  previous  one,  but  we  had  to  consider  that,  even  if  the 
output  profiles  T{y)  and  S(y)  of  the  Garwood  OBL  model  can  be  approximated  as 
piecewise  linear,  the  output  sound-speed  profile  c(y)  cannot  be  considered  as  piecewise 
linear.  This  is  because  c(y)  is  a  highly  non-linear  function  of  T.  S,  and  y.  This 
characteristic  can  be  demonstrated  by  using  a  simple  equation  like  the  one  derived  by 
Coppens  [Ref  12] 

c(T,S,y)  =   1449.05  +  4.57T  -  0.0521T-  +  0.00023T^  (4.5) 

-^  (1.333  -  0.0126T)(S  -  35)  +  0.0163y  . 

This  may  be  rewritten  as 

c  =  a  +  bT  ^  cT^  +  dT^  +  eS  +  fTS  +  gy  .  (4.6) 

Differentiation  of  c  versus  depth  y  yields 

dc  dy  =  g  4-  (b  +  2cT+3dT2  +  fS)  dT/dy  +  (e  +  f^)  dS.dy  .  (4.7) 

Hence,  in  a  particular  slab  i,  even  if  (dT  dy)j  and  (dS  dy)^  are  constant,  it  can  be  seen 
from  (4.7)  that  (dc,  dy).  is  not  a  constant  because  T  and  S  are  not  constant  in  the  slab, 
but  are  linear  functions  of  depth  y. 

Finally,  the  acoustic  model  developed  in  Chapter  II  permits  the  use  of  XBT 
profiles  discretized  with  a  depth  increment  of  one  meter  or  more.  Use  of  climatology  is 
often  unrealistic  due  to  multiple  averaging  and  too  coarse  vertical  resolution. 
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V.  APPLICATION  OF  THE  COUPLED  MODEL  TO  A  SPECIFIC  AREA 

A,       GENERALITIES 

In  order  to  investigate  the  typical  results  we  can  obtain  from  this  coupled-model, 
we  shall  simulate  some  cases  in  the  Western  Mediterranean  Sea.  The  chart  presented 
in  Figure  5.1  shows  the  site  to  be  studied:  (p  =  42.5^N,  X  =  007.5''E. 
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Figure  5.1     The  northern  region  of  the  Algero-Provencal  Basin. 

The  Mediterranean  Sea  is  one  of  the  major,  deep  closed  seas  of  the  world.  The 
water  mass  interacts  with  the  open  ocean  only  through  the  Strait  of  Gibraltar,  which 
has  a  sill  depth  (300  m)  much  shallower  than  the  average  depth  of  the  rest  of  these 
closed  waters.  The  water  below  this  sill  remains  almost  homogeneous  in  temperature 
and  salinity.    The  western  Mediterranean  is  an  area  having  very  specific  oceanographic 
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and  atmospheric  features.  Figure  5.1  depicts  the  north  part  of  the  Algero-Provencal 
Basin,  which  is  one  of  the  two  large  deep  basins  of  the  Western  Mediterranean.  The 
basin  bottom  contains  wide,  flat  abyssal  plains,  averaging  2700  meters  in  depth.  The 
water  mass  is  clearly  bounded  by  land  and  islands,  and  this  sea  is  well  ventilated  due  to 
the  quantity  of  deep  water  formed  at  the  surface  by  winter  cooling.  As  described  in 
[Rcf  13],  during  winter,  the  Alpine  relief  brings  intrusions  of  polar  continental  air 
masses  into  contact  with  the  surface  waters  of  the  northern  part  of  the  Mediterranean. 
Strong,  cold,  dry,  continental  winds  (mistral)  can  blow  for  a  few  days,  cooling  the 
surface  water  and  leading  to  a  vigorous  convective  mixing.  Such  a  strong  sinking 
motion  of  surface  water  can  mix  and  cool  uniformly  as  much  as  the  upper  thousand 
meters  of  the  water  column.  The  sea-state  is  often  characterized  by  a  short  wavelength 
swell  which  causes  a  very  rough  sea  during  storms.  Another  interesting  result  of  the 
interaction  between  the  atmosphere  and  the  sea  is  that  the  Mediterranean  water  is 
characterized  by  high  temperatures  (13'C  at  1000m)  and  high  salinities  (38.4-38.45 
g;kg)  which  are  surpassed  only  in  the  Red  Sea.  This  is  due  to  the  large  heat  input  and 
an  excess  of  evaporation  over  fresh  water  input  (precipitation  and  river  contribution). 
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Figure  5.2     Seasonal  variation  of  the  temperature  profile 
in  the  northern  part  of  the  Western  .Mediterranean. 

Figure  5.2  from  [Ref  14]  gives  an  example  of  the  seasonal  variations  observed  in  the 
Ligurian  Sea,  that  is,  the  northern  part  of  the  Algero-Provencal  Basin.  It  can  be 
observed  that,  in  this  region,  the  water  column  is  completely  isothermal  in  winter.  A 
warmer  surface  layer  evolves  progressively  during  springtime.    In  autumn,  it  decreases 
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in  temperature  hut  increases  in  depth  until  winter  conditions  are  reached  again.  The 
summer  temperature  profile  is  characterized  by  a  top  layer  about  20  to  40  meters  thick 
overlying  a  thermocline  with  large  temperature  gradients  (up  to  1  to  2'C  m  ). 

The  Mediterranean  Sea  is  also  of  interest  for  the  numerous  naval  operations 
conducted  in  that  area.  For  this  reason,  it  is  of  interest  to  study  some  of  the  features 
of  the  acoustic  propagation  that  are  related  to  atmospheric  factors.  For  example,  in 
this  region,  the  heating  by  the  sun  of  the  upper  layers  of  water,  together  with  an 
absence  of  mixing  by  the  wind,  causes  a  strong  near-surface  negative  sound-speed 
gradient  to  develop  during  the  spring  and  summer  months.  This  thermocline  overlies 
isothermal  water  at  greater  depths.  The  result  is  a  strong  shallow  internal  sound 
channel  (SOFAR  channel:  sound  fixing  and  ranging)  with  its  axial  depth  near  150 
meters.  During  the  autumn,  the  profile  returns  to  its  winter  time  conditions,  with 
isothermal  water  and  positive  sound-speed  gradient  extending  to  the  surface  of  the 
v.-ater  column,  leading  to  a  half-duct  type  of  propagation.  During  the  summer  season 
the  near-surface  negative  gradient  and  the  resulting  strong  downward  refraction  greatly 
limit  the  detection  range  of  surface  ship  hull-mounted  sonars.  By  way  of 
compensation,  the  summer  time  channel  in  the  Vlediterranean  produces  convergence 
zones  for  a  near-surface  source  similar  to  the  deep-ocean  sound  channel,  although  the 
range  of  the  zone  is  much  less  because  of  the  smaller  vertical  extent  of  the  channel 
compared  to  the  open  ocean.  Therefore,  from  the  point  of  view  of  sound  propagation, 
the  important  factor  here  is  the  existence  of  a  predominantly  isothermal  and  isohaline 
water  mass  at  depth,  which  results  in  a  sound-speed  profile  characterized  by  a  deep 
constant-gradient  section.  Only  the  upper  section  of  the  profile  will  contain 
irregularities  as  depicted  in  Figures  5.3  and  5.4.    See  [Refs.  11,14]  for  more  details. 

B.       SUMMARY  OF  THE  DIFFERENT  PARAMETERS  USED  IN  THIS  STUDY 

In  the  simulations  discussed  in  the  next  sections,  we  used  the  parameters  listed  in 
Table  3.  Four  diilerent  cases  are  to  be  analyzed,  corresponding  to  the  months  of 
December,  February',  June,  and  September.  The  boundar}'  conditions  are  assumed  to 
be  climatological  and  were  obtained  from  [Refs.  16,17,18,19].  A  precipitation  rate  of  2 
mm  hour  over  12  hours  will  be  assumed  for  hea\T  rain  cases.  Values  of  20  to  50°  o  of 
0  will  be  used  to  simulate  the  reduction  of  radiation  incident  on  the  sea  surface  due  to 
the  presence  of  clouds.  The  initial  profiles  of  temperature  and  salinity  come  from  the 
clim^atology  used  by  the  Fleet  Numerical  Oceanographic  Center  (FNOC).  Finally,  we 
will  avoid  the  complication  of  bottom  reflection  by  assuming  an  infinitely  deep  ocean. 
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Figure  5.3     SSP  for  the  Ligurian  Sea,  Western  Mediterranean 
Februarv  to  Julv  from  left  to  richt. 
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Figure  5.4     SSP  for  the  Ligurian  Sea,  Western  Mediterranean 
August  to  January  from  right  to  left. 
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C.       DECEMBER  CASE 

I.  Source  at  10  meters 

This  case  simulates  a  surface-ship  hull-mounted  sonar  and  we  will  consider 
strong  wind  (30  knots)  as  well  as  light  wind  (5  knots)  conditions. 


TABLE  3 

HEAT  FLUXES  CLIMATOLOGY  FOR  THE  MEDFrERR^ANEAN 

(W.M-) 

December 

Februar}- 

June 

September 

Qs 

130 

190 

250 

190 

Q 

400 

600 

780 

600 

Qh^Qe 

190 

120 

80 

120 

Qb^Qh  +  Qe 

280 

190 

130 

190 

Qnet 

-150 

0 

120 

0 

EV(mm;hr) 

0.3 

0.2 

0.1 

0.2 

a.   Strong  wind  conditions  (30  knots) 

The  evolution  of  the  mixed  layer  depth  h  for  different  cases  is  shown  in 
Table  4.  In  the  case  of  no  clouds  and  no  rain,  we  fmd  that  the  mixed  layer  depth 
increases  significantly  due  to  the  combined  cfTects  of  the  strong  winter  time  surface 
cooling  (Q     =-150W  nr)  plus  the  efTect  of  the  wind  stirring. 

Figure  5.5  shows  the  elTect  of  this  mixing  and  cooling  on  the  acoustic  ray 
paths.  The  first  graph  shows  the  acoustic  propagation  at  time  =  0  hour  for  a 
beamwidth  0  =  6\  and  the  second  graph  shows  the  propagation  after  24  hours  of 
strong  wind.  By  counting  the  rays,  we  fmd  that  39%  of  the  acoustic  energy  has  been 
redistributed.  More  rays  are  filling  a  deeper  mi.xed  layer,  and  the  convergence  zone 
(CZ)"'  has  weakened.  There  is  no  effect  on  the  CZ  range  which  remains  between  IS 
and  20  km. 

Figure  5.6  shows,  with  more  detail,  the  transformation  of  a  bundle  of 
refracted-surface  reflected  (R-SR)  rays  (initial  angles  between  87°  and  89.3°)  after  24 
hours  of  strong  wind.    For  a  total  of  24  R-SR  rays,  13  rays  are  now  trapped  in  the 


■^The  term  convergence  zone  (CZ)  is  appUed  to  a  phenomenon  of  focalization  of 
the  underwater  sound  in  annular  domains. 
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Figure  5.5     December,  strong  wind,  source  at  10  meters:  t  =  0,24  hrs,  6  =  6^ 
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(W)    HldlQ 


(U)    Hld^Q 


Figure  5.6    December,  transformation  of  R-SR  rays,  strong  wind:  t  =  0,24hrs. 
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mixed  layer,  increasing  the  possibility  of  detection  on  a  shallow  target  for  a  surface- 
ship  sonar.  On  the  second  graph  of  Figure  5.6,  the  weakening  of  the  CZ  is  quite 
apparent. 


TABLE  4 

EVOLUTION  OF  THE  MIXED  LAYER  DEPTH 

IN  DECEMBER  WITH 

A  STRONG  WIND 

no  rain 

no  rain 

no  rain 

heav7  rain  for  12  hrs 

no  clouds 

50%  Q^ 

20%  Q^ 

20°/o  Q^ 

t(hr) 

h(m) 

h(m) 

h(m) 

h(m) 

0 

37.0 

37.0 

37.0 

37.0 

12 

51.7 

52.5 

53.1 

50.6 

24 

60.3 

61.2 

61.7 

59.4 

36 

65.7 

67.0 

67.7 

65.7 

48 

70.8 

72.0 

72.8 

70.9 

60 

74.3 

76.0 

76.9 

75.1 

72 

78.2 

79.5 

80.4 

78.7 

From  Table  4,  we  can  see  that  adding  clouds  to  the  simulation  increases 
the  surface  cooling  and  thus  accelerates  a  somewhat  the  deepening  of  the  mixed  layer 
due  to  the  added  production  of  turbulence  by  wind  stirring.  The  eflect  of 
superimposing  heavy  rain  for  12  hours  is  a  slight  damping  of  the  mixed  layer  deepening 
rate.  The  main  conclusion  of  these  last  test  cases  is  that,  with  strong  winds,  the  overall 
efTect  of  clouds  and  rain  is  not  too  significant  and  the  general  shape  o[  the  acoustic  ray 
trace  does  not  change. 

b.    Light  wind  conditions  (5  knots) 

Light  wind  conditions  allow  us  to  analyze  the  interesting  effect  of  heavy 
rain  on  the  acoustic  propagation.  The  evolution  of  the  mixed  layer  depth  is  shown  in 
Table  5  . 

Figures  5.7  and  5.8  show  the  efiect  of  a  heavy  rain  on  a  winter  time  sound- 
speed  profile. 

The  gradient  of  the  sound-speed  versus  depth  can  be  expressed  as  : 
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Figure  5.7     December,  light  wind,  heavy  rain  for  12hrs, 
source  at  10  meters:  t  =  0,6  hrs,  6  =  5^ 
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(W)   Hld3a 


(U)    Hld3a 


Figure  5.8     December,  light  wind,  heav7  rain  for  I2hrs, 
source  at  10  meters:  t=  12,14  hrs,  6=5°. 
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dc  dy  =  cc  dl  dT  dy  +^c  ^S  dS  dy  +  dc:dp  dp;dy  (5.1) 

where  the  average  values  of  the  partial  derivatives  arc  the  following  : 

dcol  =  4.0  m  scc'°C  ,  (5.2) 


5c/t5S  =1.2  m  sec  g;kg  ,  (5.3) 

and 

5c/c'p  =  0.016  m  seem  .  (5.4) 


TABLE  5 

EVOLUTION  OF  THE  MIXED  LAYER  DEPTH  IN  DECEMBER 

WITH  LIGHT  WINDS  AND  HEAVY  RAIN  FOR  12  HRS 

t(hr) 

h(m) 

t(hr) 

h(m) 

0 

37.0 

24 

39.9 

6 

7.9 

36 

43.0 

12 

10.4 

48 

46.4 

13 

14.4 

60 

49.3 

14 

34.5 

72 

52.5 

The  temperature  term  of  the  gradient  dc/dy  is  dominant,  and  the  salinity 
term  is  usually  so  small  that  it  is  very  often  neglected.  In  the  case  of  heavy  rain, 
however,  the  term  cc  dS  dS,  dy  is  not  negligible  and  can  be  very'  important  and  even 
dominant.  In  the  mixed  layer.  dT,dy  =  0.  Thus  if  heavy  rain  occurs,  we  get 
dS,  dy>  >  0,  and  dc/dy  can  be  positive  and  much  greater  than  usual.  This  yields  a 
sharp  "inversion"  as  depicted  by  the  sound-speed  profile  of  Figures  5.7  and  5.8. 
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Another  way  to  examine  the  efTect  of  the  rain  is  to  consider  the  surface 
buoyancy  flux  Bq  : 

In  winter  time,  the  first  term  of  Bq  is  on  average  strongly  negative.  This  corresponds 
to  an  upward  buoyancy  ilux.  During  hea\7  rain  fall,  EV-PR  is  negative,  leading  to  a 
positive  second  term  in  Bq.  This  fact  can  be  verified  by  investigating  the  last  two 
columns  of  Table  4.  Rain  fall  damps  the  upward  buoyancy  flux  and,  thus,  decreases  the 
mixed  layer  deepening. 

In  Figures  5.7  and  5.S.  we  see  the  effect  of  heavw  rain  fall  on  the  acoustic 
propagation  at  time  0,  6,  12,  and  14  hours.  In  these  graphs,  the  beamwidth  has  been 
set  to  5".  After  6  hours  of  heavy  rain,  a  high  density  of  R-SR  rays  have  been  trapped 
in  the  mixed  layer.  After  12  hours,  the  efiect  is  even  much  more  important.  Almost  all 
the  rays  are  trapped  in  a  shallow  mixed  layer,  and  the  convergence  zone  has 
disappeared.  At  that  time,  the  only  chance  for  a  surface-ship  sonar  to  get  any 
detection  is  for  a  shallow  target. 

However,  this  effect  does  not  last  verv'  long  after  the  rain  stops.  We  can 
see  in  Table  5  that,  in  only  two  hours,  the  mixed  layer  depth  almost  returns  to  its 
initial  value.  We  can  also  notice  in  Figures  5.7  and  5.8  that,  only  two  hours  after  the 
end  of  the  rain,  the  acoustic  propagation  is  almost  the  same  as  for  the  initial  case  at 
t  =  0.  This  is  due  to  the  strong  winter  time  and  night  time  surface  cooling,  even  though 
the  wind  is  light.  This  fact  also  demonstrates  the  strong  effect  of  surface  cooling  on 
the  mixed  layer  deepening  in  winter  time. 
2.  Source  at  40  meters 

Now,  we  are  going  to  study  the  case  where  the  source  is  just  below  the  initial 
mixed  layer  (37m). 

a.   Strong  wind  conditions  (30  knots) 

In  the  case  of  no  clouds  and  no  rain.  Figure  5.9  shows  the  evolution  of  the 
acoustic  propagation  over  one  day  of  strong  wind  for  a  beamwidth  0  of  6°. 

Initially,  at  time  t  =  0,  all  the  rays  are  R-SR,  and  we  have  a  perfect  case  of 
convergence  zone  (CZ)  detection  with  a  strong  focusing  effect  at  the  smaller  range  of 
the  CZ  (around  15km).  After  24  hours  ol'  strong  wind,  the  deepening  efiect  of  the 
mixed    layer    redistributes    42%    of  the    acoustic    energy    into    a    mixed    layer    type 
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(W)  Mid^a 


[W)    Hlcl3a 


Figure  5.9     December,  strone  wind,  source  at  40  meters:  t  =  0,24  hrs,  0  =  6*' 
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propagation.    The  CZ  range  stays  the  same,  but  the  CZ  is  weakened  and  the  focusing 
effect  is  greatly  reduced. 

Figure  5.9  shows  that  the  detection  ability  of  a  towed  sonar  at  40  m  is 
much  improved  because  of  the  deepening  of  the  mi.xed  layer  mainly  due  to  the  wind 
stirring.  For  the  purpose  of  highlighting  this  effect,  Figure  5.10  extracts  only  the 
inOuenced  ray  bundle,  that  is,  the  bundle  with  a  beamwidth  of  2.5''.  All  the  rays  of  this 
bundle  are  transformed  from  RR  (refracted-refracted)  type  to  mixed  layer  trapped. 

If  we  would  have  added  clouds  to  the  simulation,  with  or  without  rain,  we 
would  have  found  the  same  general  shape  for  the  acoustic  ray  trace.  The  fact  confirms 
that,  in  winter  time,  the  wind  mixing  has  a  dominating  effect  on  the  acoustic 
propagation. 

b.   Light  wind  conditions  (5  knots) 

For  a  source  at  40  meters,  even  with  light  wind  conditions,  acoustic 
propagation  is  not  ver\'  sensitive  to  heavy  rainfall.  Figure  5. 1 1  shows  the  ray  trace  at 
time  0  and  after  12  hours  of  rain.  The  difference  in  shape  between  the  two  graphs  is 
not  ver\'  great.  Thus,  the  effect  of  heav}-  rain  is  to  be  considered  significant,  even  with 
light  wind  conditions  (which  is  not  usually  the  case  in  winter  time  during  storms 
passing  by),  for  only  a  shallow  transmitter,  as  seen  in  the  previous  section. 

3.  Source  at  100  meters 

If  the  source  is  sufficiently  deep  relative  to  the  mixed  layer  depth,  the 
atmospheric  forcing  has  no  effect  on  the  acoustic  propagation.  Adding  clouds  and  rain 
will  not  change  the  general  shape  of  the  acoustic  ray  paths.  Figure  5.12  depicts  the 
case  of  strong  wind  conditions  with  no  clouds  and  no  rain  at  t  =  0  and  24  hours.  For  a 
beamwidth  of  6',  the  RR  rays  stay  trapped  in  the  shallow  SOFAR  channel  and  do  not 
enter  the  mixed  layer.  This  is  because  of  the  strong  thermocline  {strong  negative 
gradient  dc  dy)  underlying  the  mixed  layer.  For  rays  to  enter  the  mixed  layer,  we 
would  need  steep  initial  angles  at  the  source,  which  is  not  realistic  for  the  usual  sonar 
beamwidth. 

In  Figure  5.13,  with  a  source  at  80  meters,  an  interesting  effect  on  the  ray 
paths  is  depicted.  The  rays  almost  appear  to  be  reflected  on  the  the  sharp 
discontinuity  created  at  the  top  of  the  thermocline  because  of  the  strong  wind  mixing. 

4.  Source  on  the  SOFAR  axis  (150  meters) 

Figure  5.14  shows  the  acoustic  ray  paths  obtained  when  the  depth  of  the 
source  is  set  to  150  meters,  that  is,  the  optimum  depth  for  the  internal  sound-channel 
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Figure  5.10     December,  strong  wind,  source  at  40  meters:  t  =  0.24  hrs,  0=2.5° 
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(W)    Hld3a 


(W)   Hlcl3Q 


Figure  5.11     December,  light  wind,  heav7  rain,  source  at  40  meters:  t  =  0,12  hrs,  6  =  6" 
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Figure  5.12     December,  strong  wind,  source  at  100  meters:  t  =  0,24  hrs,  0  =  6° 
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Figure  5.13     December,  strong  wind,  source  at  80  meters:  t  =  0,24  hrs,  Pn  =  9r  to  93* 
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Figure  5.14     December,  strong  wind,  source  at  150  m:  t  =  0,48  hrs,  6  =  9' 


69 


propagation.  Strong  wind  conditions  were  simulated  for  48  hours  and  the  beamwidth 
was  set  to  9^.  The  interesting  acoustic  propagation  is  due  to  the  sharp  positive  and 
negative  gradients  dc/dy  under  and  overlying  the  depth  of  the  SOFAR  axis.  As 
mentioned  previously,  the  source  is  deep  enough  for  the  acoustic  propagation  not  to  be 
influenced  by  the  atmospheric  forcing.  A  weak  CZ  appears  at  ranges  between  15  to  18 
km,  and  these  ranges  are  influenced  almost  not  at  all  by  the  boundar}'  conditions. 

In  that  case,  if  we  superimpose  the  second  graphs  of  Figures  5.5  and  5.14,  we 
can  see  the  combination  of  the  ray  traces  obtained  with  a  source  at  10  meters  and  a 
source  at  150  meters.  The  water  mass  is  almost  completely  filled  wiih  acoustic  energy, 
leading  to  a  high  probability  of  detection  of  a  target  whose  depth  would  be  between  0 
and  300  meters. 

5.  Source  at  500  meters 

Figure  5.15  displays  the  case  where  the  wind  has  been  blowing  for  24  hours 
with  no  clouds  and  no  rain.  We  would  have  obtained  the  same  general  shape  by 
var}'ing  the  wind  and  (or)  adding  clouds  and  rain.  As  in  the  two  previous  cases,  the 
source  is  deep  enough  for  the  acoustic  propagation  not  to  be  influenced  by  the 
atm.ospheric  factors.  In  any  case,  these  graphs  were  provided  to  show  the  excellent 
shape  of  the  acoustic  ray  trace  in  the  case  of  a  perfect  surface  reflection. 

6.  Conclusions  for  December 

On  short  time  scale,  the  atmospheric  conditions  have  an  elTect  on  the  acoustic 
propagation  only  for  a  shallow  source  (10  to  40  meters  in  our  simulations),  that  is.  for 
a  source  in  the  initial  mixed  layer  or  just  below. 

For  such  a  shallow  source,  the  wind  has  a  dominating  effect  over  rain  and 
clouds.  Mainly,  this  efiect  is  the  transformation  of  R-SR  rays  to  mixed  layer  trapped 
rays  and  the  weakening  of  the  convergence  zone. 

The  effect  of  the  surface  cooling  is  also  important  in  winter  time.  The  net 
upward  heat  flux  was  set  to  -150  W,  nr  in  our  simulations.  This  leads  to  an  average 
deepening  of  the  mixed  layer  untill  the  month  of  Februar\'  when  the  entire  water 
column  is  mixed.  This  was  depicted  in  the  light  wind  conditions  cases,  and  it  yields  a 
haif-duct  type  of  acoustic  propagation. 

The  effect  of  heavw  rain  has  been  also  depicted  for  light  wind  conditions  and 
leads  to  a  strong  shallow  inversion,  trapping  most  of  the  rays  emanating  from  a 
shallow  source  in  a  shallow  surface  duct. 
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Figure  5.15     December,  strong  wind,  source  at  500  m:  t  =  0,24  hrs,  0=6° 
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D.       FEBRUARY  CASE 

I.  Source  at  10  meters 

a.   Strong  wind  conditions  (30  knots) 

Table  6  shows  the  evolution  of  the  mixed  layer  depth  for  different  cloud 
and  rain  conditions.  In  all  the  cases,  the  mixed  layer  deepens  mainly  because  of  the 
wind  stirring  and,  as  we  saw  previously  in  the  December  case,  adding  clouds  to  the 
simulation  increases  the  deepening  of  the  mixed  layer,  while  the  rain  reduces  slightly 
this  effect. 


TABLE  6 

EVOLUTION  OF  THE  MIXED  LAYER  DEPTH 

IN  FEBRUARY  WITH  A  STRONG  WIND 

no  rain 

no  rain 

no  rain 

hea\7  rain  for  12  hrs 

no  clouds 

50%  Q^ 

20°.'o  Q^ 

20ro  Q^ 

t(hr) 

h{m) 

h(m) 

h(m) 

h(m) 

0 

70.0 

70.0 

70.0 

70.0 

12 

77.6 

79.5 

80.5 

76.8 

24 

85.5 

87.3 

88.3 

84.9 

36 

89.1 

92.3 

94.1 

90.8 

48 

94.5 

97.7 

99.5 

96.3 

60 

96.9 

101.3 

103.8 

100.7 

72 

101.2 

105.6 

108.0 

104.7 

The  initial  mixed  layer  depth  extracted  by  the  OBL  model  is  70  meters.    In 
this  model,  the  criterion  on  which  is  based  the  choice  of  this  depth  is  the  following  : 

|Abl  =  iagAT  -  PgASi  >  0.04  cm/s^  (5.6) 

where  the  symbol  A  represents  the  variation  of  buoyancy,  temperature,  and  salinity  for 
an  increment  of  one  meter  depth.  A  higher  value  for  the  criterion  (0.05  for  instance) 
would  have  increased  all  the  values  of  the  mixed  layer  depth,  but  the  evolution  would 
have  had  the  same  trend  as  the  case  depicted  in  Table  6. 
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However,  independent  of  the  choice  of  the  criterion  for  the  mixed  layer 
depth,  the  sound-speed  profile  of  Figure  5.16  shows  that  the  water  column  is  mixed 
almost  from  the  top  to  the  bottom  with  an  average  sound-speed  gradient  dc  dy  of  0.02 
sec  .  Thus,  in  February,  the  main  characteristic  of  the  acoustic  propagation  is  a  half- 
duct  type  of  propagation.  In  the  case  illustrated  by  Figure  .5.16,  the  beamwidth  has 
been  set  to  12''.  and  the  two  acoustic  ray  traces  show  the  situation  at  tmie  0  and  24 
hours.  The  wind  does  not  alTect  the  acoustic  propagation  since  the  water  mass  is 
already  almost  completly  mixed  at  t  =  0.  Nearly  identical  ray  path  shapes  would  have 
been  obtained  by  assuming  clouds  and  rain. 
b.    Lij^ht  wind  conditions  {5-10  knots) 

Table  7,  for  light  wind,  no  clouds,  and  no  rain  conditions,  highlights  the 
diurnal  variability  of  the  mixed  layer  depth  and  illustrates  how  the  surface  waters  of  the 
sea  warm  during  the  daytime  hours  on  a  sunny  day,  shallowing  the  mixed  layer  depth. 
At  night,  the  deepening  and  cooling  mixed  layer  returns  to  a  state  nearly  the  same  as 
the  initial  one.  These  diurnal  changes  have  a  profound  elTect  on  sound  transmission 
for  a  surface-ship  sonar  as  depicted  in  Figures  5.17  and  5.18. 

In  Figure  5.17,  we  simulate  a  light  wind  of  5  knots  with  no  clouds  and  no 
rain.  The  ray  traces  are  shown  at  time  instant  0  and  12  hours  for  a  beamwidth  of  8°. 
At  the  end  of  the  daytime  (t=  12  hours),  the  general  shape  of  the  acoustic  ray  trace  is 
almost  the  same  as  the  one  at  the  end  of  the  night  time  (t  =  0  hour),  except  for  the 
ranges  less  than  6  km  in  shallow  water.  This  is  depicted  in  Figure  5.18  at  t=S  hours, 
that  is.  in  the  middle  of  the  afternoon,  for  a  beamwidth  of  6".  This  reduction  of 
surface-ship  echo  ranging  on  a  shallow  target  is  a  phenomenon  often  called  the 
"afternoon  efiect." 

Next,  hea\7  rain  is  included  for  12  hours  with  light  wind  conditions.  The 
effect  is  depicted  at  time  0  and  12  hours  in  Figure  5.19.  where  the  wind  speed  was  set 
equal  10  knots  and  the  beamwidth  to  6\  At  the  initial  time,  the  acoustic  energy 
propagates  following  a  half-duct.  At  time  12  hours,  we  find  that  62'Vo  of  the  rays  are 
trapped  in  a  shallow  surface  duct  because  of  the  strong  discontinuity  in  the  sound- 
speed  profile  due  to  the  decrease  of  salinity  in  the  upper  layer  of  the  sea  caused  by  a 
weak  wind  stirring  and  damping  of  the  turbulence  by  the  strong  downward  surface 
buoyancy  Ilux  associated  with  the  rainfall.  However,  the  overall  shape  of  the  acoustic 
ray  trace  is  not  changed  significantly.  There  will  be  an  improved  possibility  of 
detection  for  a  shallow  target.  If  the  wind  continues  to  blow  at  10  knots,  the  effect 
remains  even  12  hours  after  the  rain  has  stopped. 

73 


(W)   Hld3a 


(U)    Hid3a 


Figure  5.16     Februar}',  strong  wind,  source  at  10  meters:  t  =  0,24  hrs,  6=12° 
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(U)    Hld3a 


(U)  Hid:]a 


Figure  5.17     February,  light  wind,  source  at  10  meters:  t  =  0,12  hrs,  9  =  8'^ 
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(W)    Hld3a 
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Figure  5.18     February,  light  wind,  source  at  10  meters;  t  =  0,8  hrs,  0  =  6° 
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Figure  5.19     Februar>',  light  wind,  heavy  rain,  source  at  10  meters:  t  =  0,12  hrs,  0  =  8' 
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2.  Source  at  200  meters 

In  this  case,  we  simulate  the  performance  of  a  sonar  towed  at  a  depth  of  200 
meters.  In  order  to  demonstrate  that  the  atmospheric  forcing  has  no  effect  on  sound 
transmission  from  a  deep  source.  Figure  5.20  shows  the  case  with  a  heavy  rain  for  12 
hours  together  with  a  Ught  wind  of  10  knots.  We  can  see  that  rain  does  not  influence 
the  acoustic  propagation.  Whatever  the  atmospheric  conditions  in  February,  a  deep 
towed  sonar  is  able  to  detect  any  target  moving  between  0  and  300  meters.  Thus,  the 
half-duct  type  of  propagation  is  one  of  the  best  configurations  for  underwater  detection 
by  deep  transmitters. 


TABLE  7 

EVOLUTION  OF  THE  MIXED  LAYER  DEPTH 

IN  FEBRUARY  WITH  A  LIGHT  WIND 

wind  5  knots 

wind  5  knots 

wind  10  knots 

no  rain 

heavw  rain 
for  12  hrs 

heavy  rain 
for  12  hrs 

no  clouds 

20%  Q^ 

20%  Q^ 

5 

t(hr) 

h{m) 

t(hr) 

h{m) 

t(hr) 

h(m) 

0 

70.0 

0 

70.0 

0 

70.0 

12 

7.6 

12 

2.1 

12 

12.9 

24 

69.0 

16 

3.3 

20 

14.0 

36 

7.6 

20 

7.9 

24 

16.S 

48 

69.1 

24 

67.8 

36 

33.7 

60 

7.6 

36 

70.9 

48 

70.8 

72 

69.3 

48 

73.2 

60 

73.0 

3.  Conclusions  for  February 

As  in  December,  the  atmospheric  conditions  do  not  influence  the  acoustic 
propagation  for  rays  emanating  from  a  deep  towed  sonar. 

Compared  to  December,  the  SOFAR  channel  has  disappeared,  and  there  is  a 
completely   asymmetrical,   surface-bounded  channel  (or  half-duct)  inside  which  the 
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Figure  5.20    Februar>',  light  wind,  heavy  rain,  source  at  200  meters:  t  =  0,12  hrs,  6  =  8' 
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propagation  is  made  of  R-SR  rays  only.  Because  of  this  half-duct  type  of  propagation, 
February-  is  one  of  the  best  months  in  the  Mediterranean  for  underwater  detection  by 
both  deep  and  shallow  transmitters.  This  generalisation  not  take  in  consideration  the 
sea  state  which  can  generate  noise  and  reflection  loss  when  the  sea  is  rough,  because 
our  acoustic  model  assumes  a  perfect  surface  reflection. 

With  light  wind  conditions,  surface-ship  echo  ranging  can  be  poorer  in  the 
afternoon  on  a  shallow  target  because  of  the  "afternoon  effect"  due  to  surface  water 
heating. 
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Figure  5.21     Variation  of  sound-channel  width  d  for  similar  temperature 
changes  in  (a)  the  Mediterranean  and  (b)  the  Atlantic. 


E.       JUNE  CASE 

l-'igures  5.3  and  5.4  show  that  the  shape  of  the  sound-speed  profile  stays 
approximatively  the  same  during  the  entire  summer  season,  that  is,  from  June  to 
September.  Therefore,  in  this  section,  we  will  just  simulate  the  atmospheric  conditions 
and  the  acoustic  propagation  in  June.  All  of  our  observations  will  be  applicable  to  the 
other  months  of  summer  time  even  though  the  atmospheric  factors  change  somewhat, 
as  illustrated  in  Table  3. 

In  summer,  so  much  heat  is  added  to  the  upper  layer  of  the  Mediterranean  Sea 
that  negative  gradients  in  temperature  and,  thus,  in  sound-speed  develop  in  the  upper 


SO 


layer  of  the  sea.  This  is  underlain  by  isothermal  water  at  a  relatively  shallow  depth. 
The  consequence  is  that  a  sound  channel  exists  in  summer  in  the  Mediterranean  and  is 
characterized  as  follows: 

The  sound  channel  is,  in  all  cases,  strongly  asymmetrical.  The  sound-speed 
profile  contains  a  sharp  discontinuity  between  sound-speed  gradients  of  completely 
ditTerem  orders  of  miagnitude:  0.017  sec  on  the  lower  side,  due  to  the  pressure  elTect, 
and  -1  to  -5  sec      on  the  upper  side  in  the  thermocline. 


TABLE  8 

EVOLUTION  OF  THE  MIXED  LAYER  DEPTH 
IN  JUNE  WITH  STRONG  WINDS 

no  rain 

no  rain 

no  rain 

heavy  rain  for  12  hrs 

no  clouds 

50%  Q^ 

20%  Q^ 

20%  Q^ 

t{hr) 

h(m) 

h{m) 

hmi) 

h(m) 

0 

2.0 

2.0 

2.0 

2.0 

12 

20.3 

20.8 

21.1 

20.3 

24 

25.7 

26.2 

26.5 

25.8 

36 

29.4 

30.3 

30.8 

30.0 

48 

32.9 

33.9 

34.5 

33.6 

60 

35.5 

36.9 

37.7 

36.8 

72 

1 

38.3 

39.8 

40.7 

39.8 

The  sound  channel  axis  is  ver\'  close  to  the  surface,  and  the  curvature  of  the 
sound-speed  profile  near  this  minimum  is  verv'  great. 

The  width  of  the  sound  channel  is  smaller  than  in  the  Atlantic  and  is  extremely 
variable  with  the  temperature  of  the  surface  v.-ater,  as  depicted  in  Figure  5.21.  The 
usual  width  of  the  sound  channel  in  summer  is  on  the  order  of  1400  to  1800  meters. 
The  depth  of  the  Mediterranean  is  sufiicient  for  this  channel  not  to  be  intercepted  by 
the  sea  bottom  in  most  cases. 
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1.  Source  at  10  meters 

a.   Strong  wind  conditions  (30  knots) 

The  evolution  of  the  mixed  layer  depth  h  for  different  cases  is  shown  in 
Table  8.  In  the  case  of  no  clouds  and  no  rain,  we  find  that  the  mixed  layer  depth 
increases  significantly  during  the  first  12  hours  because  of  the  strong  wind  mixing  in 
spite  of  the  strong  summer  time  heating.  As  seen  in  the  previous  sections,  for  short 
time  scales,  the  wind  is  the  doininating  factor  influencing  both  the  mixed  layer  depth 
and  the  acoustic  propagation. 

Figure  5.22  shows  the  acoustic  ray  trace  at  time  0  for  a  beamwidth  8=8°. 
After  24  hours  of  30  knot  winds,  the  mixed  layer  deepens  to  25.7  meters.  Even  if  the 
convergence  zone  weakens  slightly,  the  range  of  the  CZ  is  not  affected.  However,  the 
acoustic  propagation  changes  greatly  for  the  less  steep  rays.  Thirty  percent  of  the  rays 
are  now  trapped  in  a  shallow  mixed  layer,  improving  slightly  the  probability  of 
detection  using  a  surface-ship  sonar.  However,  the  unfavorable  conditions  for 
underwater  detection  for  a  shallow  source  is  well  depicted  in  Figure  5.23. 


TABLE  9 

MAXIMUM  AND  MINIMUM  ML  DEPTHS  AND  TEMPERATURES 
IN  JUNE  FOR  LIGHT  WIND  CONDITIONS  WITH  NO  CLOUDS 

day  2 

9a.m. 

h  =  4.8m 

T=20.13"C 

5  p.m. 

h=l.lm 

T=21.3UC 

day  3 

9a.m. 

h=5.4m 

T=20.15^C 

5p.m. 

h=  1.1m 

T=21.33'C 

day  4 

9  a.m. 

h=  5.6m 

T=20.21°C 

5p.m. 

h=l.lm 

T=21.39°C 

b.    Light  wind  conditions  (5  knots) 

Table  9  shows  the  evolution  of  the  maximum  and  minimum  mixed  layer 
depths  and  high  temperatures  rise  during  light  wind  conditions  with  no  clouds.  The 
mixed  layer  responds  to  the  diurnal  cycle  of  solar  radiation,  shallowing  in  daytime 
because  of  solar  heating,  and  deepening  at  night  because  of  surface  cooling.    We  notice 
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(W)    Hld3a 


(u)  Hidia 


Figure  5.22    June,  strong  wind,  source  at  10m:  t  =  0,24  hrs,  6=8' 
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Figure  5.23    June,  strong  wind,  source  at  10m:  t  =  0,24  hrs,  0=8' 
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that  the  general  trend  of  the  minin^.um  and  maximum  temperature  is  to  increase, 
leading  progressively  to  the  extreme  sound-speed  profile  of  August  (see  Figure  5.4). 
Xethertheless,  the  mixed  layer  stays  sufTiciently  shallow  not  to  influence  the  acoustic 
propagation,  even  for  a  shallow  transnutter. 

Adding  clouds  to  the  simulation  (20'^o  Q,)  only  deepens  the  mixed  layer  to 
9  meters  after  48  hours,  which  has  almost  no  effect  on  the  acoustic  propagation.  A 
longer  period  of  cloud  covering  would  not  be  a  realistic  simulation  for  summertime  in 
the  VIediterranean. 

Finally,  heavy  rain  keeps  the  mixed  layer  shallow  (around  1.5  meters),  and 
sound  propagation  is  not  influenced. 

2.  Source  on  the  SO  FAR  a.xis  (90  meters) 

Figure  5.24  simulates  the  case  where  the  source  is  on  the  shallow  sound- 
channel  axis  characterizing  the  Mediterranean  summer  time.  The  beamwidth  is  set  to 
8°.  The  strong  asymmetry  of  the  sound  channel  is  readily  apparent.  As  previously 
mentioned,  the  wind  does  not  influence  the  sound  channel  propagation  for  6  =  8".  The 
rays  at  the  source  are  not  sufTiciently  steep  to  be  affected  and,  even  if  the  sound 
channel  width  decreases,  the  propagation  does  not  change. 

However,  in  order  to  study  the  effect  due  to  the  decrease  of  the  sound  channel 
width.  Figure  5.25  shows  the  ray  traces  for  steep  initial  angles  P=8r,  82%  83',  97°, 
98",  and  99"  at  time  0  and  after  48  hours  of  strong  winds  (30  knots).  The  channel 
width  decreases  from  1400  to  900  meters.  The  steepest  rays  (p  =  8r,  82",  98".  and  99") 
are  transformed  from  RR  type  to  R-SR  type,  but  the  influence  of  this  phenomenon  on 
the  acoustic  propagation  is  almost  not  perceptible  with  an  acoustic  ray  tracing  model 
as  shown  by  Figure  5.25.  A  transmission  loss  model  would  be  necessary'  to  account  for 
the  transmission  loss  due  to  the  reflection  on  a  rough  sea  surface. 

3.  Conclusions  for  June 

In  June,  and  for  the  reminder  of  the  summer,  the  sound  propagation  for  a 
deep  transmitter  is  characterized  by  a  sound-channel  type  of  propagation  and  is  not 
influenced  by  the  atmospheric  forcing  on  a  short  time  scale  of  a  few  days. 

Because  of  the  absence  of  a  mixed  layer  and  the  presence  of  a  sharp 
thermocline  due  to  the  important  effect  of  surface  heating,  the  detection  capabilities  of 
a  surface-ship  sonar  are  poor  and  limited  to  the  convergence  zone  with  ranges  between 
30  and  35  km.  However,  strong  winds  can  create  a  shallow  mixed  layer  in  which 
shallow  targets  could  be  detected  as  depicted  in  Figure  5.26,  where  rays  emanating 
from  both  shallow  and  deep  transmitters  are  traced. 
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Figure  5.24    June,  strong  wind,  source  at  90  meters:  t  =  0,24  hrs,  G  =  8' 
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(W)    Hld3a 


(w)  ^lAia 


Figure  5.25    June,  strong  wind,  source  at  90  meters,  t  =  0,48  hrs,  P  =  81^82^83^97^98^99' 
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(U)   Hld3a 


(W)    Hld3a 


Figure  5.26 


June,  strong  wind,  sources  at  10  and  90  meters 
t  =  0,48  hrs,  e=8°. 
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VI.  CONCLUSION 

On  a  short  time  scale  of  a  few  hours  to  three  days,  the  atmospheric  conditions 
influence  the  acoustic  propagation  only  for  a  shallow  source.  The  ray  trace  emanating 
from  a  deep  transmitter  is  almost  not  affected  by  the  atmospheric  forcing.  This 
conclusion  is  independent  of  the  type  of  propagation,  depending  on  the  time  of  the 
year  (deep-sound  channel,  convergence  zone,  or  half  duct). 

The  wind  is  the  dominant  atmospheric  factor.  For  short  time  scales,  the  coupled 
model  is  not  very  scnsitiN'e  to  the  other  boundary  conditions  (solar  irradiation,  and 
rainfall).  This  is  a  favorable  finding  in  the  sense  that  the  wind  is  the  easiest  parameter 
to  measure  and  to  forecast  at  sea.  However,  the  importance  of  a  heavy  rainfall  has 
been  demonstrated  in  the  case  of  light  wind  conditions. 

The  initial  temperature  profile  plays  a  determinant  role  and  has  to  be  provided  as 
accurately  as  possible.  The  coupled  model  is  not  very  sensitive  to  the  initial  salinity 
profile. 

A.  ADVANTAGES  OF  THE  COUPLED  MODEL 

The  graphic  output  sequence  from  the  coupled  model  shows  qualitatively  the 
acoustic  propagation  and  its  evolution  under  the  forcing  of  specified  atmospheric 
conditions,  and  an  assumed  initial  temperature  profile. 

The  simplicity  of  the  model  makes  it  possible  to  be  implemented  on  a  simple 
desktop  computer  with  some  graphics  capability. 

A  current  XBT,  digitized  by  using  a  one-meter  depth  increment,  can  be  used  as 
input  in  the  model. 

B.  WEAKNESS  OF  THE  COUPLED  MODEL 

The  fact  that  the  coupled  model  does  not  account  ver\'  well  for  the  operating 
frequency  is  inherent  to  a  ray  tracing  routine.  The  user  of  such  a  model  will  have  to 
simulate  a  beamwidth  at  the  source  corresponding  approximatively  to  a  gi\en 
frequency. 

The  output  of  the  coupled  model  does  not  provide  such  quantitative  results  as 
median  detection  range,  convergence  zone  range,  S,  N  ratio  at  the  receiver,  and  it  does 
not  account  for  scattering  and  attenuation. 
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Perfect  surface  reflection  has  been  assumed  whatever  the  sea-state.  An 
infmitively-deep  ocean  has  been  assumed  in  order  to  avoid  the  complication  of  bottom 
reflection.  However,  the  model  can  handle  perfect  bottom  reflection  or  absorption. 
Furthermore,  the  earth  curvature  has  been  neglected. 

The  coupled  model  assumes  horizontal  homogeneity  of  the  water  mass,  which  we 
considered  to  be  realistic  for  the  horizontal  ranges  considered. 

Finally,  we  neglected  horizontal  advection,  which  is  compatible  with  the  short 
time  scales  studied,  especially  in  the  Mediterranean  where  the  currents  are  weak. 

C.       RECOMMENDATIONS 

Further  simulations  could  be  applied  to  different  locations  in  the  Atlantic, 
Pacific,  or  Indian  Oceans,  where  climatology  is  more  available. 

Real  data  could  be  used  to  initialize  the  coupled  model  and  the  boundary' 
conditions  could  be  computed  from  observed  atmospheric  and  oceanic  conditions.  The 
resulting  outputs  of  the  model  could  be  compared  with  corresponding  simulations. 

Then,  one  might  couple  the  OBL  model  with  one  of  the  different  acoustic  models 
used  by  the  US  Navy  (such  as  RAYMODE,  FACT,  or  PE  models)  in  order  to  obtain 
quantitative  results  including  MDR  or  CZR  for  usual  operating  frequencies  and  a 
variety  of  environmental  conditions. 
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APPENDIX 
RAY  TRACING  SUBROUTINE  USE  EXAMPLE 

File  giving  an  example  of  sound  speed  profile  to  be  used  with  RAY 


0. 

1482. 

30. 

14S2.9 

55. 

1484.1 

90. 

1480.5 

150. 

1482. 

225. 

1484.1 

300. 

1486.2 

LCDR  JACQUES  FOURNIOL.  FRENCH  NAVY. 
PARAMETER  (Y0=50. ,B0=90. , DB= . 5 ,M=7 ,MM=6 ,N=300 ,NN=5000 ,RANGE=30 . ) 
DIMENSION  BET(M) ,CC(0:MM) ,YY(0:MM) 
DIMENSION  C(0:N) ,G(-1:N) ,YC(0:N) 
DIMENSION  Y(0:NN) ,Z(0:NN) 
READ(05,10)(YY(I) , CC(I ) , 1=0 ,MM) 
10    FORMAT(F7.2,2X,F7.2) 

CALL  RAY(YO,BET,BO,DB,M,CC,YY,MM,C,G,YC,N,Y,Z,NN,RANGE) 

STOP 

END 


C 

c 
c 
c 
c 
c 

C 

c 
c 
c 

C 

c 
c 
c 
c 
c 
c 
c 
c 
c 
c 
c 


SUBROUTINE  USING  THE  MM+1  POINTS  OF  A  STRAIGTH  LINE  SEGMENT 
CONTINUOUS  PROFILE  (  CC  VS  YY  )  TO  TRACE  M  ACOUSTIC  RAYS  ISSUED 
FROM  A  SOURCE  AT  A  DEPTH  YO . 

REF.  "UNDERWATER  ACOUSTICS-A  LINEAR  SYSTEMS  THEORY  APPROACH" 
BY  ZIOMEK  L.  P237-23S.  ACADEMIC  PRESS,  ORLANDO,  FLORIDA.  1985. 
PERFECT  SURFACE  AND  BOTTOM  REFLECTIONS.  FLAT  BOTTOM. 


BETAO 

M 

BET 

BO 

DB 

YO 

MM+1 

CC 

YY 

NN 

Y 

Z 

N 

C 


STARTING  ANGLE  OF  RAY,  REF  VERTICAL,  IN  RAD. 

NB  OF  RAYS  SENT. 

ARRAY  OF  THE  INITIAL  ANGLES  IN  DEGREES  AT  THE  SOURCE, 

INITIAL  ANGLE  OF  THE  UPPER  RAY  OF  THE  BUNDLE. 

INCREMENT  OF  INITIAL  ANGLES  IN  BET. 

SOURCE  DEPTH  IN  METER. 

NB  OF  PTS  IN  PROVIDED  SSP. 

PROVIDED  SSP  IN  M/SEC. 

DEPTH  OF  THE  PTS  OF  THIS  SSP  IN  METERS. 

INDEX  OF  RANGE. 

DEPTH  IN  METERS  AND 

HOR.  DISTANCE  TRAVELLED  BY  A  RAY  IN  KM. 

VAX   INTEGER  VALUE  OF  DEPTH  IN  METER. 

SSP  .1  METER  INCREMENT. 
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C        G      :  GRADIENT  OF  SSP. 

C        YC     :  ARRAY  OF  DEPTH  FOR  PLOTTING  SSP. 

C        RANGE  :  MAX.  RANGE  OF  THE  PLOTIN  KM. 

C        K      :  DEPTH  INDEX. 

C     N  MUST  BE  THE  INTEGER  VALUE  OF  THE  MAX.  DEPTH  YY(MM)  OF  THE 

C     INITIAL  SSP. 

C 

SUEROUT INE  RAY ( YO , BET , BO , DB , M , CC, YY , MM , C, G , YC, N , Y , Z , NN , RANGE ) 
C 

DIMENSION  BET(M) ,CC(0:MM) ,YY(0:MM) 
DIMENSION  C(0:N) ,G(-1 :N) ,YC(0:N) 
DIMENSION  Y(0:NN) ,Z(0:NN) 
DATA  PI/3.141592654/ 
Y(0)=YY(0) 
YC(0)=YY(0) 
JMIN=0 
C(0)=CC(0) 
C 

C     COMPUTE  THE  SSP  C  AND  ITS  GRADIENT  G  WITH  1  METER  DEPTH  INCREMENT 
DO  10  1=1, MM 

DYY=YY(I)-YY(I-1) 
JMAX=IFIX(YY(I)) 
GG=(CC(I)-CC(I-1))/DYY 
DO  30  J=JMIN,JMAX-1 
G(J)=GG 
C(J+1)=C(J)+GG 
30       YC(J+1)=FL0AT(J+1) 

JMIN=JMAX 
10    CONTINUE 
C 

C     CREATE  GRADIENT  AT  SFC  AND  BOTTOM  TO  HANDLE  REFLECTION 
G(-1)=G(0) 
G(N)=G(N-1) 
C 

C      DO  03  J=1,M 
C03    BET(J)=90.-FLOAT(J)*0.5 
DO  03  J=1,M 
03    BET(J)=B0-FLOAT(J)*DEL3 
C 

C     RESEARCH  OF  CMIN  AND  CMAX  FOR  PLOTTING 
CMIN=1550. 
CMAX=1450. 
DO  40  1=0, N 
CMIN=AMIN1 { CMIN , C ( I ) ) 
40    CMAX=AMAX1 ( CMAX , C ( I ) ) 
CMIN=AINT(CMIN) 
CMAX=AINT(CMAX)+1. 
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c 
c 

C     PLOT  THE  SSP  AND  PREPARE  PLOT  FOR  RAYS 
C      CALL  COMPRS 

CALL  0X41(4107) 
C      CALL  TEK618 

C      CALL  SHERPA( 'THESEOUT' , 'B' ,1) 
C       CALL  HWROT ( ' AUTO ' ) 

CALL  BLOWUP (0.745) 

CALL  PAGE (11. ,5.5) 

CALL  PHYSOR( 1.2,1.) 

CALL  AREA2D(1. ,3.) 

CALL  THKFRM(.Ol) 

CALL  XNAME('SSP  (M/S) $ ' , -100) 

CALL  YNAMEC  DEPTH  (M)$',100) 

CALL  CROSS 

CALL  YAXAMG(0.) 

CALL  YTICKS(5) 

CALL  XTICKS(IFIX(CMAX-CMIN)) 

CALL  GRAF(CMIN,CMAX-CMIN,CMAX,FLOAT(N) ,-50. ,0. ) 

CALL  CURVE ( C , YC , M+ 1 , 0 ) 

CALL  FRAME 

CALL  ENDGR(O) 

CALL  PHYS0R(2.7,1.) 

CALL  AREA2D(7.3,3.) 

CALL   XNAME  (  '  RAI^IGE    (KM)  $',-100) 

CALL  YNAME( '$',100) 

CALL  YNONUM 

CALL  CROSS 

CALL  XTICKS(5) 

CALL  GRAF(0. ,5. , RANGE , FLOAT (N) ,-50. ,0.) 
C 

DO  05  J=1,M 

Y(0)=Y0 

Z(0)=0. 

K=IFIX(Y0) 

BETA0=BET(J)/180.'^PI 

CO=C(K) 

B=SIN(BETAO)/CO 
C 

c 

DO  20  1=1, NN 

NNN=I 
C 
C     CHECK  IF  SFC  REFLECTION 

IF  (K.EQ.O)  BETAO=PI-BETAO 
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CHECK  IF  BOTTOM  REFLECTION 
IF  (K.EQ.N)  BETAO=PI-BETAO 

IF  ((G(K).GT.O).AND.(G(K-l).GT.O.))  THEN 
IF  (BETAO.GT.PI/2.)  THEN 
K=K-1 

BETA=PI-ASIN(B'^C(K)) 
Y(I)=FLOAT(K) 

Z ( I ) =Z ( I - 1 ) + ( COS (BETAO ) - COS ( BETA) ) / (B*G (K) ) 
BETAO=BETA 
ELSE 

BETAC=ASIN(C(K)/C(K+1)) 
IF  ( BETAO. LT.BETAC)  THEN 
K=K+1 

BETA=ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z( I )=Z(I-1)+(C0S (BETAO) -COS (BETA) )/(B*G(K-l)) 
BETA0=3ETA 
ELSE 

Y(I)=FLOAT(K) 

Z(I)=Z(I-1)+2.'^C0S  (BETAO)/ (B*G(K)) 
BETAO=PI-BETAO 
END  IF 
ENDIF 
ELSE 

IF  ((G(K).LT.O.).AND.(G(K-l).LT.O.))  THEN 
IF  (BETAO. LT. PI/2.)  THEN 
K=K+1 

EETA=ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z( I )=Z(I-1)+(C0S (BETAO) -COS (BETA) )/(B*G(K-l)) 
BETAO=BETA 
ELSE 

BETAC=PI-ASIN(C(K)/C(K-1)) 
IF  ( BETAO. GT.BETAC)  THEN 
K=K-1 

BETA=PI-ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z(I)=Z(I-1)+(C0S (BETAO) -COS (BETA) )/(B*G(K)) 
BETAO=BETA 
ELSE 

Y(I)=FLOAT(K) 

Z(I)=Z(I-1)+2.*COS(BETAO)/(B*G(K-1)) 
BETAO=PI -BETAO 
ENDIF 
ENDIF 
ELSE 
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IF((G(K).LT,0.) .AND. (G(K-1 ) .GT . 0 . ) )  THEN 
IF  (BETAO.GT.PI/2.)  THEN 
K=K-1 

BETA=PI-ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z(I)=Z(I-1)+(COS(BETAO)-COS(BETA))/(B*G(K)) 
BETAO=BETA 
ELSE 
K=K+1 

BETA=ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z(I)=Z(I-1)+(COS(BETAO) -COS (BETA) )/(B*G(K-l)) 
BETAO=BETA 
END  IF 
ELSE 

IF  ((G(K) .GT.O.) .AND.(G(K-l).LT.O.))  THEN 
IF  (EETAO.GT.PI/2.)  THEN 

BETAC=PI-ASIN(C(K)/C(K-1)) 
IF  (BETAO.GT.BETAC)  THEN 
K=K-1 

BETA=PI-ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z(I)=Z(I-1)  +  (COS(BETAO)-COS(BETA))/(B'^G(K)) 
BETAO=BETA 
ELSE 

Y(I)=FLOAT(K) 

Z(I)=Z(I-l)+2.*COS(BETA0)/(B*G(K-l)) 
BETAO=PI-BETAO 
END  IF 
ELSE 

BETAC=ASIN(C(K)/C(K+1)) 
IF  (BETAO.LT.BETAC)  THEN 
K=K+1 

BETA=ASIN(B*C(K)) 
Y(I)=FLOAT(K) 

Z(I)=Z(I-1)+(COS(BETAO)-COS(BETA))/(B*G(K-1) 
BETAO=BETA 
ELSE 

Y(I)=FLOAT(K) 

Z(I)=Z(I-l)+2.*COS(BETA0)/(B*G(K)) 
BETAO=PI-BETAO 
END  IF 
ENDIF 
END  IF 
ENDIF 
ENDIF 
ENDIF 
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IF(Z(I).GT.RANGE*1000.)  GO  TO  21 

20  CONTINUE 
C 

c 

21  CONTINUE 
C 

C      CONVERSION  M  TO  KM 
DO  50  1=0, NNN 
50    Z(I)=Z(I)*0.001 
C 

c 

C     PLOT  THE  RAYS 

CALL  CURVE (Z,Y,NNN+1,0) 
C 

05    CONTINUE 
C 

CALL  FRAME 
C      CALL  HEADIN( 'ACOUSTIC  RAY  TRACING$ ',100 , 2 . , 1 ) 

CALL  ENDPL(O) 

CALL  DONEPL 

STOP 

END 
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